esMethods_ManuscriptAnalysis
Karen Sasmita
11/27/2020
Setup
Load necessary packages and clear workspace
library(ggplot2)
library(viridis)
library(lme4)
library(lmerTest)
library(emmeans)
library(knitr)
library(dplyr)
library(drc)
library(aomisc)
library(nlme)
library("ggpubr")
library(broom)
library(gridExtra)
library(gtable)
library(grid)
library(psycho)
rm(list = ls())Setup plot themes
#Set text size and colors
txt.size = 12
txt.color = "black"
txt.face = "bold"
###define specs for agreement values plots over sample sizes
#agreement value colors
d.color = "#FFC700" #actual data
random.color = "#00B885" #random data
crossMov.color = "#4D4D4D" #agreement values when actual data compared to cross movie.
factor.order <- c('crossMov', 'random', 'sample')
factor.color <- c(crossMov.color, random.color, d.color)
factor.shape <- c(24, 23, 21) #triangle, diamond, circle
factor.label <- c('Cross-movie', 'Random', 'Same')
#jitter for individual agreement values
jitter.alpha = .1
jitter.width = 0.2
#point for mean agreement value
point.size = 5
point.alpha = 1
#errorbar width
errorbar.width = .8
errorbar.alpha = 1
#elbow related plot
elbow.size = 8
elbowPos.scale = .01
#define specs for curvature function fittings.
vir_7 = viridis(n = 7) #load the first 8 viridis colors
line.colors <- c(vir_7[1],vir_7[1], vir_7[2],vir_7[3], vir_7[4], vir_7[5], vir_7[6]) #Linear, Exponential, logarithmic, Asymptotic, Power Curve, Yield loss.
dot.color <- vir_7[7] #color for dotplot
fit.jitter.alpha = .05
fit.jitter.width = 0.3
fit.line.size = 1.5
fit.line.alpha = 1
###define specs for saving agreement and function fit plots
figure.width = 25
figure.height = 19
figure.unit = "cm"
figure.bg = "transparent"
figure.filetype = "png"
grain.labs <- c("Coarse", "Fine")
names(grain.labs) <- c("c", "f")
sampleSize = seq(from = 2, to = 32, by = 2)
#Set theme for plotting
theme.esMethods <- theme(
panel.grid.minor = element_blank(),
panel.grid.major = element_blank(),
#panel.background= element_blank(),
panel.border = element_blank(),
panel.spacing = unit(.05,'in'),
panel.background = element_rect(fill = "transparent",colour = NA),
plot.background = element_rect(fill = "transparent",colour = NA),
axis.line = element_line(size = .2, colour = "black"),
axis.title = element_text(size = txt.size),
axis.text = element_text(size = txt.size),
strip.background = element_rect(colour = NA, fill = "transparent"),
strip.text = element_text(size = txt.size),
legend.title = element_blank(),
legend.key = element_rect(colour = NA, fill = NA),
plot.title = element_text(hjust = 0.5)
)Define necessary functions
Function: elbow_finder Adapted from https://stackoverflow.com/questions/2018178/finding-the-best-trade-off-point-on-a-curve/2022348#2022348
Function description: to find the stabilization point (elbow) of change in segmentation agreement based on best fitting growth/ decay functions.
Returns: coordinates of curve that is furthest from the straight line connecting the minimum and maximum points.
Parameter x_values: array of x coordinates of curve
Precondition x_values: numeric
Parameter x_values: array of y coordnates of curve
Precondition y_values: numeric
elbow_finder <- function(x_values, y_values) {
# Max values to create line
max_x_x <- max(x_values)
max_x_y <- y_values[which.max(x_values)]
min_x_x <- min(x_values)
min_x_y <- y_values[which.min(x_values)]
max_df <- data.frame(x = c(min_x_x, max_x_x), y = c(min_x_y, max_x_y))
# Creating straight line between the max values
fit <- lm(max_df$y ~ max_df$x)
# Distance from point to line
distances <- c()
for(i in 1:length(x_values)) {
distances <- c(distances, abs(coef(fit)[2]*x_values[i] - y_values[i] + coef(fit)[1]) / sqrt(coef(fit)[2]^2 + 1^2))
}
# Max distance point
x_max_dist <- x_values[which.max(distances)]
y_max_dist <- y_values[which.max(distances)]
return(c(x_max_dist, y_max_dist))
}Function: functionFits adapted from https://www.statforbiology.com/nonlinearregression/usefulequations
Function description: fit raw agreement values to growth (for peakiness, agreement index, surprise index, detection accuracy and predictive value) or decay (for peak-to-peak distance) functions
Parameter y: array of dependent variable
Precondition y: numeric
Parameter x: array of independent variable
Precondition x: factor (note: this parameter is a factor because it was taken from data frame that’s fed into lmer)
Parameter sampSize: sequence of sample size
Precondition sampSize:
functionFit <- function(y, x, sampSize){
x <- as.numeric(as.character(x))
linear <- lm(y~x)
exponential <- tryCatch(drm(y~x, fct = DRC.expoDecay()), error = function(e) drm(y~x, fct = EXD.2()))
log <- drm(y ~ x, fct = DRC.logCurve())
asympReg <- drm(y ~ x, fct = DRC.asymReg())
powerCurve <- tryCatch(drm(y~x, fct = DRC.powerCurve()), error = function(e) drm(y~x, fct = L.4()))
yl <- nls(y~NLS.YL(x, a, i)) #yielded erratic behavior for normativeHit
models <- c('linear','exponential', 'log', 'asympReg', 'powerCurve', 'yl')
model.prediction <- data.frame()
for (model in models){
assign("prediction", get(model))
model.prediction <- rbind(model.prediction, data.frame(func = model, sampleSize = x, prediction = predict(prediction)))
}
bic <- BIC(linear,exponential,log,asympReg,powerCurve, yl)
aic <- AIC(linear,exponential, log, asympReg, powerCurve, yl)
functionFit.summary <- data.frame(AIC_df = aic$df,
AIC = aic$AIC,
BIC_df = bic$df,
BIC = bic$BIC)
row.names(functionFit.summary) <- models
mod.predicted <- predict(get(models[match(min(bic$BIC), bic$BIC)]), data.frame(sampSize))#, interval = 'confidence', level = .99
elbow_pos <- elbow_finder(sampSize, mod.predicted)
model <- models[match(min(bic$BIC), bic$BIC)]
return(list(model.prediction, functionFit.summary,elbow_pos, model, coef(get(model))[3]))
}Function: plotIndivAgreement
Function description: function to plot agreement values
Parameter df: data frame containing values for plotting
Precondition df: data frame
Parameter var: variable name for y values (i.e. which agreement measure)
Precondition var: string
Parameter ylab: y-axis label
Precondition ylab: string
Parameter ymin: minimum y axis value
Precondition ymin: numeric
Parameter ymax: maximum y axis value
Precondition ymax: numeric
Parameter title: plot title
Precondition title: string
plotIndivAgreement <- function(df, var, ylab, ymin, ymax, title){
factors <- factor.order[(length(unique(df$testCond))%%3):3]
factor.colors <- factor.color[(length(unique(df$testCond))%%3):3]
factor.shapes <- factor.shape[(length(unique(df$testCond))%%3):3]
factor.labels <- factor.label[(length(unique(df$testCond))%%3):3]
transform(df, testCond=factor(testCond,levels=factors))
elbow = mean(df$elbow)
figure <- ggplot(df, aes(x = sampleSize, y = get(var, df), fill = testCond, color = testCond, shape = testCond))+
geom_jitter(alpha = jitter.alpha, width = jitter.width)+
stat_summary(geom = "errorbar", fun.data = "mean_cl_normal", width = errorbar.width, color = 'black', alpha = errorbar.alpha)+
geom_point(stat = "summary", fun = "mean", size = point.size, alpha = point.alpha, color = 'black')+
scale_fill_manual(values = factor.colors, labels = factor.labels)+
scale_color_manual(values = factor.colors, labels = factor.labels)+
scale_shape_manual(values = factor.shapes, labels = factor.labels)+
#scale_x_continuous(breaks = seq(from = 2, to = 32, by = 2))+
geom_text(color = "black", x = elbow/2, y = ymin + (elbowPos.scale*(ymax - (ymin))), label = '^', size = elbow.size, show.legend = FALSE)+
#facet_wrap(~grain, scales='fixed', labeller = labeller(grain = c("c" = "Coarse", "f" = "Fine")))+
# geom_text(color = 'black', x = 6, y = -.2, label = "+", size = txt.size, show.legend = FALSE)+
labs(y = ylab, x = "Sample Size", title = title)+
ylim(ymin,ymax)+
theme.esMethods
#assign(paste('figure.', g, sep = ''), figure)`
#ggsave(paste(paste('figure.', g, sep = ''), ".png", sep = ""))
return(figure)
}Function: plotFunctionFits
Function description: function to plot fitted growth/ decay functions
Parameter df: data frame containing values for plotting (i.e. function fit predictions)
Precondition df: data frame
Parameter original_df: data frame containing raw agremement values (to add as jitter)
Precondition df: data frame
Parameter var: original_df variable name for agreement values (i.e. which agreement measure)
Precondition var: string
Parameter ylab: y-axis label
Precondition ylab: string
Parameter ymin: minimum y axis value
Precondition ymin: numeric
Parameter ymax: maximum y axis value
plotFunctionFits <- function(df, original_df, var, ylab, ymin, ymax){
figure <- ggplot(df,aes(x = sampleSize, y = prediction, color = func))+
geom_jitter(data = original_df[original_df$testCond == "sample",],
aes(x = as.numeric(as.character(sampleSize)), y = get(var, original_df)), alpha = fit.jitter.alpha, width = fit.jitter.width, fill = dot.color, color = dot.color)+
geom_line(size = fit.line.size, alpha = fit.line.alpha)+
scale_color_manual(values = line.colors,
labels = c("Linear", "Exponential", "Logarithmic", "Asymptotic", "Power Curve", "Yield Loss"))+
geom_point(data = original_df[original_df$testCond == "sample",], aes (x = as.numeric(as.character(sampleSize)), y = get(var, original_df)), stat = "summary", fun = "mean", shape = 21, size = 4, color = 'black', fill = dot.color)+
scale_x_continuous(breaks = c(seq(from =2, to = 32, by = 2)))+
ylim(ymin,ymax)+
labs(y = ylab, x = "Sample size")+
#facet_wrap(~grain, scales='fixed', labeller = labeller(grain = c("c" = "Coarse", "f" = "Fine")))+
theme.esMethods
return(figure)
}Analyse data
Peakiness
Load dataset and remove outliers for evAct peakiness (value > 20*sd)
setwd('../data/bootstrapped')
peakiness.comm <- read.delim("esMethods_Peakiness_adjust_c0.1_f0.05_100Iterations_commercial.txt", head = TRUE)
peakiness.comm$sampleSize <- as.factor(peakiness.comm$sampleSize)
peakiness.comm$peakiness <- peakiness.comm$act.peakiness/peakiness.comm$min.peakiness
peakiness.comm$logPeakiness <- log(peakiness.comm$peakiness, base = 10)
peakiness.evAct <- read.delim("esMethods_Peakiness_adjust_c0.1_f0.05_100Iterations_evAct.txt", head = TRUE)
peakiness.evAct$peakiness <- peakiness.evAct$act.peakiness/peakiness.evAct$min.peakiness
peakiness.evAct <- peakiness.evAct[peakiness.evAct$peakiness < (mean(peakiness.evAct$peakiness) + (20*sd(peakiness.evAct$peakiness))),]
peakiness.evAct$sampleSize <- as.factor(peakiness.evAct$sampleSize)
peakiness.evAct$logActPeakiness <- log(peakiness.evAct$act.peakiness, base = 10)
peakiness.evAct$logPeakiness <- log(peakiness.evAct$peakiness, base = 10)Fit peakiness values to growth and decay functions
##commercial
peakiness.fit.coarse.comm <- functionFit(peakiness.comm$peakiness[peakiness.comm$grain == "c" & peakiness.comm$testCond == 'sample'], peakiness.comm$sampleSize[peakiness.comm$grain == "c" & peakiness.comm$testCond == 'sample'], sampleSize)
peakiness.fit.fine.comm <- functionFit(peakiness.comm$peakiness[peakiness.comm$grain == "f" & peakiness.comm$testCond == 'sample'], peakiness.comm$sampleSize[peakiness.comm$grain == "f" & peakiness.comm$testCond == 'sample'], sampleSize)
peakiness.functionFit <- peakiness.fit.coarse.comm[[1]] %>% mutate(grain = "c")
peakiness.functionFit <- rbind(peakiness.functionFit, peakiness.fit.fine.comm[[1]] %>% mutate(grain = "f"))
peakiness.functionFit$grain <- as.factor(peakiness.functionFit$grain)
#peakiness.functionFit.comm.figure <- plotFunctionFits(peakiness.functionFit, peakiness.comm[peakiness.comm$testCond == "sample",], "peakiness", "Peakiness (raw)")
#add elbow value to main dataframe
for (i in 1:nrow(peakiness.comm)){
if(peakiness.comm$grain[i] == 'c'){
peakiness.comm$elbow[i] = 0 #coarse commercial peakiness did not stabilize
} else {
peakiness.comm$elbow[i] = peakiness.functionFit[[3]][1]
}
}
##everyday Activities
peakiness.fit.coarse.evAct <- functionFit(peakiness.evAct$peakiness[peakiness.evAct$grain == "c" & peakiness.evAct$testCond == "sample"], peakiness.evAct$sampleSize[peakiness.evAct$grain == 'c' & peakiness.evAct$testCond == "sample"], sampleSize)
peakiness.fit.fine.evAct <- functionFit(peakiness.evAct$peakiness[peakiness.evAct$grain == "f" & peakiness.evAct$testCond == "sample"], peakiness.evAct$sampleSize[peakiness.evAct$grain == 'f' & peakiness.evAct$testCond == "sample"], sampleSize)
peakiness.evAct.functionFit <- peakiness.fit.coarse.evAct[[1]] %>% mutate(grain = "c")
peakiness.evAct.functionFit <- rbind(peakiness.evAct.functionFit, peakiness.fit.fine.evAct[[1]] %>% mutate(grain = "f"))
peakiness.evAct.functionFit$grain <- as.factor(peakiness.evAct.functionFit$grain)
#peakiness.functionFit.evAct.figure <- plotFunctionFits(peakiness.evAct.functionFit, peakiness.evAct[peakiness.evAct$testCond == "sample",], "peakiness", "Peakiness (raw)") + ylim(0,50)
#add elbow value to main dataframe
for (i in 1:nrow(peakiness.evAct)){
if(peakiness.evAct$grain[i] == 'c'){
peakiness.evAct$elbow[i] = peakiness.fit.coarse.evAct[[3]][1]
} else {
peakiness.evAct$elbow[i] = 0 #fine evAct peakiness did not stabilize
}
}Plot peakiness over sample size
peakiness.comm.coarse.plot <- plotIndivAgreement(peakiness.comm[peakiness.comm$grain == 'c',], "logPeakiness", "Peakiness (log)", 0, 2.5, "Coarse")
peakiness.comm.fine.plot <- plotIndivAgreement(peakiness.comm[peakiness.comm$grain == 'f',],"logPeakiness", "Peakiness (log)", 0, 2.5, "Fine")
peakiness.evAct.coarse.plot <- plotIndivAgreement(peakiness.evAct[peakiness.evAct$grain == 'c', ], "logPeakiness", "Peakiness (log)", 0, 2.5, "")
peakiness.evAct.fine.plot <- plotIndivAgreement(peakiness.evAct[peakiness.evAct$grain == 'f', ], "logPeakiness", "Peakiness (log)", 0, 2.5, "")
peakiness.plot.combined <- ggarrange(peakiness.comm.coarse.plot, peakiness.comm.fine.plot,
peakiness.evAct.coarse.plot, peakiness.evAct.fine.plot,
nrow = 2, ncol = 2,
labels = c("A","", "B", ""),
common.legend = TRUE, legend = "bottom",
align = "hv")
peakiness.plot.combined
#ggsave("esMethod_peakiness.pdf", device = cairo_pdf, width = 25, height = 19, units = "cm", bg = "transparent")
#ggsave("esMethod_peakiness.eps", device = cairo_ps, width = 25, height = 19, units = "cm", bg = "transparent")Plot peakiness function fits
peakiness.comm.coarse.funcPlot <- plotFunctionFits(peakiness.functionFit[peakiness.functionFit$grain == 'c',], peakiness.comm[peakiness.comm$grain == 'c' & peakiness.comm$testCond == 'sample',], "peakiness", "Peakiness", 2, 24)
peakiness.comm.fine.funcPlot <- plotFunctionFits(peakiness.functionFit[peakiness.functionFit$grain == 'f',], peakiness.comm[peakiness.comm$grain == 'f' & peakiness.comm$testCond == 'sample',], "peakiness", "Peakiness", 2, 24)
peakiness.evAct.coarse.funcPlot <- plotFunctionFits(peakiness.evAct.functionFit[peakiness.evAct.functionFit$grain == 'c',], peakiness.evAct[peakiness.evAct$grain == 'c' & peakiness.evAct$testCond == 'sample',], "peakiness", "Peakiness", 2, 50)
peakiness.evAct.fine.funcPlot <- plotFunctionFits(peakiness.evAct.functionFit[peakiness.evAct.functionFit$grain == 'f',], peakiness.evAct[peakiness.evAct$grain == 'f' & peakiness.evAct$testCond == 'sample',], "peakiness", "Peakiness", 2, 50)
peakiness.funcFit.plot.combined <- ggarrange(peakiness.comm.coarse.funcPlot, peakiness.comm.fine.funcPlot,
peakiness.evAct.coarse.funcPlot, peakiness.evAct.fine.funcPlot,
nrow = 2, ncol = 2,
labels = c("A","", "B", ""),
common.legend = TRUE, legend = "bottom",
align = "hv")
peakiness.funcFit.plot.combined
#ggsave("esMethod_peakiness_functionFit.pdf", device = cairo_pdf, width = 25, height = 19, units = "cm", bg = "transparent")Peakiness linear models and analysis
peakiness.comm.lmer <- lmer(peakiness~sampleSize*grain*testCond + (1|movieName), data = peakiness.comm)
#emm_options(pbkrtest.limit = 3000)
peakiness.comm.samp_vs_rand <- summary(emmeans(peakiness.comm.lmer, c("testCond"), by = "grain", contr = list(c1 = c(-1, 1))), infer = T)
peakiness.comm.pairwise <- summary(emmeans(peakiness.comm.lmer, "testCond",by = c("sampleSize", "grain"), weights = "proportional", contr = "revpairwise", adjust = "none")$contrasts, infer = TRUE)#random vs sample pairwise/ sample size
peakiness.comm.poly <- summary(contrast(emmeans(peakiness.comm.lmer,'sampleSize', by = c( 'grain', 'testCond')), 'poly', max.degree = 3), infer = TRUE)#Polynomial contrast
peakiness.evAct.lmer <- lmer(peakiness~sampleSize*grain*testCond + (1|movieName), data = peakiness.evAct)
#emm_options(lmerTest.limit = 25600)
peakiness.evAct.samp_vs_rand <- summary(emmeans(peakiness.evAct.lmer, c("testCond"), by = "grain", contr = list(c1 = c(-1, 1))), infer = T)# overall random vs sample
peakiness.evAct.pairwise <- summary(emmeans(peakiness.evAct.lmer, "testCond",by = c("sampleSize", "grain"), weights = "proportional", contr = "revpairwise", adjust = "none")$contrasts, infer = TRUE)#pairwise comparison
peakiness.evAct.poly <- summary(contrast(emmeans(peakiness.evAct.lmer,'sampleSize', by = c( 'grain', 'testCond')), 'poly', max.degree = 3), infer = TRUE)#Polynomial contrastPeakiness Sample - random contrast commercial-lab
peakiness.comm.samp_vs_rand$contrasts## grain = c:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio p.value
## c1 4.587464 0.03646846 Inf 4.515987 4.658941 125.793 <.0001
##
## grain = f:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio p.value
## c1 5.859044 0.03646846 Inf 5.787567 5.930521 160.661 <.0001
##
## Results are averaged over the levels of: sampleSize
## Confidence level used: 0.95Peakiness sample-random contrast everyday-online
peakiness.evAct.samp_vs_rand$contrasts## grain = c:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio p.value
## c1 15.657130 0.2507069 Inf 15.165753 16.148506 62.452 <.0001
##
## grain = f:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio p.value
## c1 6.619238 0.2506576 Inf 6.127958 7.110518 26.407 <.0001
##
## Results are averaged over the levels of: sampleSize
## Confidence level used: 0.95Peakiness smallest pairwise effect commercial-lab
peakiness.comm.pairwise %>% dplyr::filter(p.value < .05, z.ratio > 0) %>% dplyr::group_by(grain) %>% filter(z.ratio == min(z.ratio))## # A tibble: 2 x 10
## # Groups: grain [2]
## contrast sampleSize grain estimate SE df asymp.LCL asymp.UCL
## <fct> <fct> <fct> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 sample … 6 c 1.34 0.146 Inf 1.05 1.63
## 2 sample … 4 f 0.347 0.146 Inf 0.0614 0.633
## # … with 2 more variables: z.ratio <dbl>, p.value <dbl>Peakiness smallest pairwise effect everyday-online
peakiness.evAct.pairwise %>% dplyr::filter(p.value < .05, z.ratio > 0) %>% dplyr::group_by(grain) %>% filter(z.ratio == min(z.ratio))## # A tibble: 2 x 10
## # Groups: grain [2]
## contrast sampleSize grain estimate SE df asymp.LCL asymp.UCL
## <fct> <fct> <fct> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 sample … 2 c 10.7 1.00 Inf 8.75 12.7
## 2 sample … 2 f 3.68 1.00 Inf 1.72 5.65
## # … with 2 more variables: z.ratio <dbl>, p.value <dbl>Peakiness polynomial fit commercial-lab
peakiness.comm.poly[peakiness.comm.poly$testCond == "sample",]## contrast grain testCond estimate SE df asymp.LCL
## 7 linear c sample 309.806809 3.803925 Inf 302.35125
## 8 quadratic c sample 33.310526 7.795726 Inf 18.03118
## 9 cubic c sample 4.709911 103.547837 Inf -198.24012
## 10 linear f sample 336.377790 3.803925 Inf 328.92223
## 11 quadratic f sample -76.594747 7.795726 Inf -91.87409
## 12 cubic f sample -302.832594 103.547837 Inf -505.78262
## asymp.UCL z.ratio p.value
## 7 317.26236 81.44398802 0.000000e+00
## 8 48.58987 4.27292157 1.929283e-05
## 9 207.65994 0.04548536 9.637204e-01
## 10 343.83335 88.42913694 0.000000e+00
## 11 -61.31540 -9.82522292 8.767983e-23
## 12 -99.88256 -2.92456707 3.449357e-03Peakiness polynomial fit commercial-lab
peakiness.evAct.poly[peakiness.evAct.poly$testCond == "sample",]## contrast grain testCond estimate SE df asymp.LCL asymp.UCL
## 7 linear c sample 197.32892 26.16873 Inf 146.03915 248.61868
## 8 quadratic c sample -146.86548 53.63168 Inf -251.98164 -41.74933
## 9 cubic c sample 1707.65861 712.29000 Inf 311.59586 3103.72135
## 10 linear f sample -80.43273 26.14541 Inf -131.67678 -29.18868
## 11 quadratic f sample 32.04620 53.58214 Inf -72.97286 137.06526
## 12 cubic f sample 752.69029 711.71232 Inf -642.24022 2147.62081
## z.ratio p.value
## 7 7.5406381 4.676775e-14
## 8 -2.7384094 6.173718e-03
## 9 2.3974204 1.651097e-02
## 10 -3.0763620 2.095433e-03
## 11 0.5980762 5.497891e-01
## 12 1.0575766 2.902485e-01Peak-to-peak distance
Load dataset
setwd('../data/bootstrapped')
peaktopeak.comm <- read.delim("esMethods_PeakPeakDistance_c0.1_f0.05_100Iterations_commercial.txt", head = TRUE)
peaktopeak.comm$sampleSize <- as.factor(peaktopeak.comm$sampleSize)
peaktopeak.comm$logPeaktopeak <- log(peaktopeak.comm$nearestDistance, base = 10)
peaktopeak.evAct <- read.delim("esMethods_PeakPeakDistance_c0.1_f0.05_100Iterations_evAct.txt", head = TRUE)
peaktopeak.evAct$sampleSize <- as.factor(peaktopeak.evAct$sampleSize)
peaktopeak.evAct$logPeaktopeak <- log(peaktopeak.evAct$nearestDistance, base = 10)Fit peak-to-peak values to growth and decay functions
##commercial
peaktopeak.fit.coarse.comm <- functionFit(peaktopeak.comm$nearestDistance[peaktopeak.comm$grain == "c" & peaktopeak.comm$testCond == 'sample'], peaktopeak.comm$sampleSize[peaktopeak.comm$grain == "c" & peaktopeak.comm$testCond == 'sample'], sampleSize)
peaktopeak.fit.fine.comm <- functionFit(peaktopeak.comm$nearestDistance[peaktopeak.comm$grain == "f" & peaktopeak.comm$testCond == 'sample'], peaktopeak.comm$sampleSize[peaktopeak.comm$grain == "f" & peaktopeak.comm$testCond == 'sample'], sampleSize)
peaktopeak.functionFit <- peaktopeak.fit.coarse.comm[[1]] %>% mutate(grain = "c")
peaktopeak.functionFit <- rbind(peaktopeak.functionFit, peaktopeak.fit.fine.comm[[1]] %>% mutate(grain = "f"))
peaktopeak.functionFit$grain <- as.factor(peaktopeak.functionFit$grain)
#peaktopeak.functionFit.comm.figure <- plotFunctionFits(peaktopeak.functionFit, peaktopeak.comm[peaktopeak.comm$testCond == "sample",], "nearestDistance", "Peak-to-peak distance (s)")
for(i in 1:nrow(peaktopeak.comm)){
if(peaktopeak.comm$grain[i] == 'c'){
peaktopeak.comm$elbow[i] = peaktopeak.fit.coarse.comm[[3]][1]
}else{
peaktopeak.comm$elbow[i] = peaktopeak.fit.fine.comm[[3]][1]
}
}
##everyday Activities
peaktopeak.fit.coarse.evAct <- functionFit(peaktopeak.evAct$nearestDistance[peaktopeak.evAct$grain == "c" & peaktopeak.evAct$testCond == "sample"], peaktopeak.evAct$sampleSize[peaktopeak.evAct$grain == 'c' & peaktopeak.evAct$testCond == "sample"], sampleSize)
peaktopeak.fit.fine.evAct <- functionFit(peaktopeak.evAct$nearestDistance[peaktopeak.evAct$grain == "f" & peaktopeak.evAct$testCond == "sample"], peaktopeak.evAct$sampleSize[peaktopeak.evAct$grain == 'f' & peaktopeak.evAct$testCond == "sample"], sampleSize)
peaktopeak.evAct.functionFit <- peaktopeak.fit.coarse.evAct[[1]] %>% mutate(grain = "c")
peaktopeak.evAct.functionFit <- rbind(peaktopeak.evAct.functionFit, peaktopeak.fit.fine.evAct[[1]] %>% mutate(grain = "f"))
peaktopeak.evAct.functionFit$grain <- as.factor(peaktopeak.evAct.functionFit$grain)
#peaktopeak.functionFit.evAct.figure <- plotFunctionFits(peaktopeak.evAct.functionFit, peaktopeak.evAct[peaktopeak.evAct$testCond == "sample",], "nearestDistance", "Peak-to-peak distance (s)") + ylim(0,50)
for(i in 1:nrow(peaktopeak.evAct)){
if(peaktopeak.evAct$grain[i] == 'c'){
peaktopeak.evAct$elbow[i] = peaktopeak.fit.coarse.evAct[[3]][1]
}else{
peaktopeak.evAct$elbow[i] = peaktopeak.fit.fine.evAct[[3]][1]
}
}Plot peak-to-peak over sample size
Log10 transformed peak-to-peak distance over increasing sample sizes for segmentation of: A) commercial-lab and B) everyday-online.
peaktopeak.comm.coarse.plot <- plotIndivAgreement(peaktopeak.comm[peaktopeak.comm$grain == 'c',], "logPeaktopeak", "Peak to peak distance (log)", 0, 2, "Coarse")
peaktopeak.comm.fine.plot <- plotIndivAgreement(peaktopeak.comm[peaktopeak.comm$grain == 'f',], "logPeaktopeak", "Peak to peak distance (log)", -0.5, 1.5, "Fine")
peaktopeak.evAct.coarse.plot <- plotIndivAgreement(peaktopeak.evAct[peaktopeak.evAct$grain == 'c', ], "logPeaktopeak", "Peak to peak distance (log)", 0, 2, "")
peaktopeak.evAct.fine.plot <- plotIndivAgreement(peaktopeak.evAct[peaktopeak.evAct$grain == 'f', ], "logPeaktopeak", "Peak to peak distance (log)", -0.5, 1.5, "")
peaktopeak.plot.combined <- ggarrange(peaktopeak.comm.coarse.plot, peaktopeak.comm.fine.plot,
peaktopeak.evAct.coarse.plot, peaktopeak.evAct.fine.plot,
nrow = 2, ncol = 2,
labels = c("A","", "B", ""),
common.legend = TRUE, legend = "bottom",
align = "hv")
peaktopeak.plot.combined
setwd('../plots')
ggsave("esMethod_peaktopeak.pdf", device = cairo_pdf, width = 25, height = 19, units = "cm", bg = "transparent", dpi = 300)
# ggsave(plot = peaktopeak.plot.combined, file = "esMethod_peaktopeak.png",
# type = "cairo-png", bg = "transparent",
# width = 25, height = 19, units = "cm", dpi = 300)Plot peak-to-peak function fits
Peak-to-peak values fit to various growth and decay functions over increasing sample sizes for segmentation of: A) commercial-lab and B) everyday-online.
#df, original_df, var, ylab, ymin, ymax, asymptote
peaktopeak.comm.coarse.funcPlot <- plotFunctionFits(peaktopeak.functionFit[peaktopeak.functionFit$grain == 'c',], peaktopeak.comm[peaktopeak.comm$grain == 'c' & peaktopeak.comm$testCond == 'sample',], "nearestDistance", "Peak-to-peak distance (s)", 0, 15)
peaktopeak.comm.fine.funcPlot <- plotFunctionFits(peaktopeak.functionFit[peaktopeak.functionFit$grain == 'f',], peaktopeak.comm[peaktopeak.comm$grain == 'f' & peaktopeak.comm$testCond == 'sample',], "nearestDistance", "Peak-to-peak distance (s)", 0, 10)
peaktopeak.evAct.coarse.funcPlot <- plotFunctionFits(peaktopeak.evAct.functionFit[peaktopeak.evAct.functionFit$grain == 'c',], peaktopeak.evAct[peaktopeak.evAct$grain == 'c' & peaktopeak.evAct$testCond == 'sample',], "nearestDistance", "Peak-to-peak distance (s)", 0, 25)
peaktopeak.evAct.fine.funcPlot <- plotFunctionFits(peaktopeak.evAct.functionFit[peaktopeak.evAct.functionFit$grain == 'f',], peaktopeak.evAct[peaktopeak.evAct$grain == 'f' & peaktopeak.evAct$testCond == 'sample',], "nearestDistance", "Peak-to-peak distance (s)", 0, 5)
peaktopeak.funcFit.plot.combined <- ggarrange(peaktopeak.comm.coarse.funcPlot, peaktopeak.comm.fine.funcPlot,
peaktopeak.evAct.coarse.funcPlot, peaktopeak.evAct.fine.funcPlot,
nrow = 2, ncol = 2,
labels = c("A","", "B", ""),
common.legend = TRUE, legend = "bottom",
align = "hv")
peaktopeak.funcFit.plot.combined
setwd('../plots')
ggsave("esMethod_peaktopeak_functionFit.pdf", device = cairo_pdf, width = 25, height = 19, units = "cm", bg = "transparent")Peak-to-peak linear models and analysis
#commercial-lab
peaktopeak.comm.lmer <- lmer(nearestDistance~sampleSize*grain*testCond + (1|movieName), data = peaktopeak.comm)
#emm_options(pbkrtest.limit = 19200)
peaktopeak.comm.poly <- summary(contrast(emmeans(peaktopeak.comm.lmer,'sampleSize', by = c( 'grain', 'testCond')), 'poly', max.degree = 3), infer = TRUE)#Polynomial contrast
peaktopeak.comm.rand_vs_samp <- summary(emmeans(peaktopeak.comm.lmer, c("testCond"), by = "grain", contr = list(c1 = c(0, 1, -1))), infer = T)# overall random - sample
peaktopeak.comm.crossMov_vs_samp <- summary(emmeans(peaktopeak.comm.lmer, c("testCond"), by = "grain", contr = list(c1 = c(1, 0, -1))), infer = T)# overall crossMov - sample
peaktopeak.comm.crossMov_vs_rand <- summary(emmeans(peaktopeak.comm.lmer, c("testCond"), by = "grain", contr = list(c1 = c(1, -1, 0))), infer = T)# overall crossMov - random
peaktopeak.comm.pairwise <- summary(emmeans(peaktopeak.comm.lmer, "testCond",by = c("sampleSize", "grain"), weights = "proportional", contr = "pairwise", adjust = "none")$contrasts, infer = TRUE)#random vs sample pairwise/ sample size
#everyday-online
peaktopeak.evAct.lmer <- lmer(nearestDistance~sampleSize*grain*testCond + (1|movieName), data = peaktopeak.evAct)
#emm_options(pbkrtest.limit = 3000)
peaktopeak.evAct.poly <- summary(contrast(emmeans(peaktopeak.evAct.lmer,'sampleSize', by = c( 'grain', 'testCond')), 'poly', max.degree = 3), infer = TRUE)#Polynomial contrast
peaktopeak.evAct.rand_vs_samp <- summary(emmeans(peaktopeak.evAct.lmer, c("testCond"), by = "grain", contr = list(c1 = c(0, 1, -1))), infer = T)# overall random - sample
peaktopeak.evAct.crossMov_vs_samp <- summary(emmeans(peaktopeak.evAct.lmer, c("testCond"), by = "grain", contr = list(c1 = c(1, 0, -1))), infer = T)# overall crossMov - sample
peaktopeak.evAct.crossMov_vs_rand <- summary(emmeans(peaktopeak.evAct.lmer, c("testCond"), by = "grain", contr = list(c1 = c(1, -1, 0))), infer = T)# overall crossMov - random
peaktopeak.evAct.pairwise <- summary(emmeans(peaktopeak.evAct.lmer, "testCond",by = c("sampleSize", "grain"), weights = "proportional", contr = "pairwise", adjust = "none")$contrasts, infer = TRUE)#random vs sample pairwise/ sample sizePeak-to-peak random - sample contrast commercial-lab
peaktopeak.comm.rand_vs_samp$contrasts## grain = c:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio p.value
## c1 3.0491790 0.03624087 Inf 2.9781482 3.1202098 84.136 <.0001
##
## grain = f:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio p.value
## c1 0.7547491 0.03624087 Inf 0.6837183 0.8257799 20.826 <.0001
##
## Results are averaged over the levels of: sampleSize
## Confidence level used: 0.95Peak-to-peak random - sample contrast everyday-online
peaktopeak.evAct.rand_vs_samp$contrasts## grain = c:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio p.value
## c1 2.6898060 0.1011654 Inf 2.4915256 2.8880865 26.588 <.0001
##
## grain = f:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio p.value
## c1 0.3761267 0.1011654 Inf 0.1778462 0.5744071 3.718 0.0002
##
## Results are averaged over the levels of: sampleSize
## Confidence level used: 0.95Peak-to-peak cross movie - sample contrast commercial-lab
peaktopeak.comm.crossMov_vs_samp$contrasts## grain = c:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio p.value
## c1 3.3981889 0.03624087 Inf 3.327158 3.4692197 93.767 <.0001
##
## grain = f:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio p.value
## c1 0.7305018 0.03624087 Inf 0.659471 0.8015326 20.157 <.0001
##
## Results are averaged over the levels of: sampleSize
## Confidence level used: 0.95Peak-to-peak cross movie - sample contrast everyday-online
peaktopeak.evAct.crossMov_vs_samp$contrasts## grain = c:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio p.value
## c1 5.3137203 0.1011654 Inf 5.1154399 5.5120008 52.525 <.0001
##
## grain = f:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio p.value
## c1 0.7301477 0.1011654 Inf 0.5318673 0.9284282 7.217 <.0001
##
## Results are averaged over the levels of: sampleSize
## Confidence level used: 0.95Peak-to-peak cross movie - random contrast commercial-lab
peaktopeak.comm.crossMov_vs_rand$contrasts## grain = c:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio
## c1 0.34900988 0.03624087 Inf 0.27797908 0.42004068 9.630
## p.value
## <.0001
##
## grain = f:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio
## c1 -0.02424736 0.03624087 Inf -0.09527815 0.04678344 -0.669
## p.value
## 0.5035
##
## Results are averaged over the levels of: sampleSize
## Confidence level used: 0.95Peak-to-peak cross movie - random contrast everyday-online
peaktopeak.evAct.crossMov_vs_rand$contrasts## grain = c:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio p.value
## c1 2.6239143 0.1011654 Inf 2.4256338 2.8221948 25.937 <.0001
##
## grain = f:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio p.value
## c1 0.3540211 0.1011654 Inf 0.1557406 0.5523015 3.499 0.0005
##
## Results are averaged over the levels of: sampleSize
## Confidence level used: 0.95Peak-to-peak smallest sample size for significant pairwise effect commercial-lab
peaktopeak.comm.pairwise %>% dplyr::filter(p.value < .05, z.ratio > 0) %>% dplyr::group_by(grain, contrast) %>% filter(sampleSize == min(as.numeric(as.character(sampleSize))))## # A tibble: 5 x 10
## # Groups: grain, contrast [5]
## contrast sampleSize grain estimate SE df asymp.LCL asymp.UCL
## <fct> <fct> <fct> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 crossMo… 2 c 0.295 0.145 Inf 0.0111 0.579
## 2 crossMo… 2 c 2.70 0.145 Inf 2.42 2.98
## 3 random … 2 c 2.40 0.145 Inf 2.12 2.69
## 4 crossMo… 2 f 0.557 0.145 Inf 0.273 0.841
## 5 random … 2 f 0.726 0.145 Inf 0.441 1.01
## # … with 2 more variables: z.ratio <dbl>, p.value <dbl>Peak-to-peak smallest sample size for significant pairwise effect everyday-online
peaktopeak.evAct.pairwise %>% dplyr::filter(p.value < .05, z.ratio > 0) %>% dplyr::group_by(grain, contrast) %>% filter(sampleSize == min(as.numeric(as.character(sampleSize))))## # A tibble: 3 x 10
## # Groups: grain, contrast [3]
## contrast sampleSize grain estimate SE df asymp.LCL asymp.UCL
## <fct> <fct> <fct> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 crossMo… 2 c 1.25 0.405 Inf 0.457 2.04
## 2 crossMo… 2 c 3.79 0.405 Inf 3.00 4.58
## 3 random … 2 c 2.54 0.405 Inf 1.75 3.33
## # … with 2 more variables: z.ratio <dbl>, p.value <dbl>Peak-to-peak smallest pairwise effect commercial-lab
peaktopeak.comm.pairwise %>% dplyr::filter(p.value < .05, z.ratio > 0) %>% dplyr::group_by(grain, contrast) %>% filter(z.ratio == min(z.ratio))## # A tibble: 5 x 10
## # Groups: grain, contrast [5]
## contrast sampleSize grain estimate SE df asymp.LCL asymp.UCL
## <fct> <fct> <fct> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 crossMo… 2 c 0.295 0.145 Inf 0.0111 0.579
## 2 crossMo… 2 c 2.70 0.145 Inf 2.42 2.98
## 3 random … 4 c 2.10 0.145 Inf 1.82 2.38
## 4 crossMo… 2 f 0.557 0.145 Inf 0.273 0.841
## 5 random … 12 f 0.642 0.145 Inf 0.358 0.927
## # … with 2 more variables: z.ratio <dbl>, p.value <dbl>Peak-to-peak smallest pairwise effect everyday-online coarse
peaktopeak.evAct.pairwise %>% dplyr::filter(p.value < .05, z.ratio > 0) %>% dplyr::group_by(grain, contrast) %>% filter(z.ratio == min(z.ratio))## # A tibble: 3 x 10
## # Groups: grain, contrast [3]
## contrast sampleSize grain estimate SE df asymp.LCL asymp.UCL
## <fct> <fct> <fct> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 crossMo… 2 c 1.25 0.405 Inf 0.457 2.04
## 2 crossMo… 2 c 3.79 0.405 Inf 3.00 4.58
## 3 random … 4 c 2.09 0.405 Inf 1.30 2.88
## # … with 2 more variables: z.ratio <dbl>, p.value <dbl>Peak-to-peak largest pairwise effect everyday-online fine (note: non significant effect)
peaktopeak.evAct.pairwise %>% dplyr::filter(p.value > .05, z.ratio > 0) %>% dplyr::group_by(grain, contrast) %>% filter(z.ratio == max(z.ratio))## # A tibble: 3 x 10
## # Groups: grain, contrast [3]
## contrast sampleSize grain estimate SE df asymp.LCL asymp.UCL
## <fct> <fct> <fct> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 random … 22 f 0.433 0.405 Inf -3.60e-1 1.23
## 2 crossMo… 32 f 0.408 0.405 Inf -3.85e-1 1.20
## 3 crossMo… 32 f 0.793 0.405 Inf -5.30e-5 1.59
## # … with 2 more variables: z.ratio <dbl>, p.value <dbl>Peak-to-peak polynomial fit commercial-lab
peaktopeak.comm.poly[peaktopeak.comm.poly$testCond == "sample",]## contrast grain testCond estimate SE df asymp.LCL
## 13 linear c sample -215.37582 3.780185 Inf -222.78484
## 14 quadratic c sample 295.65601 7.747074 Inf 280.47203
## 15 cubic c sample -2278.86403 102.901608 Inf -2480.54747
## 16 linear f sample -66.42806 3.780185 Inf -73.83709
## 17 quadratic f sample 89.66520 7.747074 Inf 74.48122
## 18 cubic f sample -711.27868 102.901608 Inf -912.96212
## asymp.UCL z.ratio p.value
## 13 -207.96679 -56.974941 0.000000e+00
## 14 310.84000 38.163572 0.000000e+00
## 15 -2077.18058 -22.146049 1.138796e-108
## 16 -59.01904 -17.572702 3.987230e-69
## 17 104.84919 11.574073 5.577059e-31
## 18 -509.59523 -6.912221 4.771230e-12Peak-to-peak polynomial fit everyday-online
peaktopeak.evAct.poly[peaktopeak.evAct.poly$testCond == "sample",]## contrast grain testCond estimate SE df asymp.LCL
## 13 linear c sample -144.34337 10.55228 Inf -165.02546
## 14 quadratic c sample 200.56105 21.62574 Inf 158.17538
## 15 cubic c sample -2036.14458 287.24691 Inf -2599.13818
## 16 linear f sample -32.64406 10.55228 Inf -53.32615
## 17 quadratic f sample 53.53322 21.62574 Inf 11.14756
## 18 cubic f sample -451.80064 287.24691 Inf -1014.79424
## asymp.UCL z.ratio p.value
## 13 -123.66128 -13.678881 1.357655e-42
## 14 242.94671 9.274183 1.789848e-20
## 15 -1473.15098 -7.088482 1.355909e-12
## 16 -11.96197 -3.093556 1.977736e-03
## 17 95.91889 2.475440 1.330719e-02
## 18 111.19296 -1.572865 1.157501e-01Agreement Index
Load dataset
setwd('../data/bootstrapped')
agreement.comm <- read.delim("esMethods_AgreementIndex_100Iterations_commercial.txt", head = TRUE)
agreement.comm$sampleSize <- as.factor(agreement.comm$sampleSize)
agreement.evAct <- read.delim("esMethods_AgreementIndex_100Iterations_evAct.txt", head = TRUE)
agreement.evAct$sampleSize <- as.factor(agreement.evAct$sampleSize)Fit Agreement Index values to growth and decay functions
##commercial-lab
agreement.fit.coarse.comm <- functionFit(agreement.comm$agreementIndex[agreement.comm$grain == "c" & agreement.comm$testCond == 'sample'], agreement.comm$sampleSize[agreement.comm$grain == "c" & agreement.comm$testCond == 'sample'], sampleSize)
agreement.fit.fine.comm <- functionFit(agreement.comm$agreementIndex[agreement.comm$grain == "f" & agreement.comm$testCond == 'sample'], agreement.comm$sampleSize[agreement.comm$grain == "f" & agreement.comm$testCond == 'sample'], sampleSize)
agreement.functionFit <- agreement.fit.coarse.comm[[1]] %>% mutate(grain = "c")
agreement.functionFit <- rbind(agreement.functionFit, agreement.fit.fine.comm[[1]] %>% mutate(grain = "f"))
agreement.functionFit$grain <- as.factor(agreement.functionFit$grain)
for(i in 1:nrow(agreement.comm)){
if(agreement.comm$grain[i] == 'c'){
agreement.comm$elbow[i] = agreement.fit.coarse.comm[[3]][1]
} else {
agreement.comm$elbow[i] = agreement.fit.fine.comm[[3]][1]
}
}
##everyday Activities-lab
agreement.fit.coarse.evAct <- functionFit(agreement.evAct$agreementIndex[agreement.evAct$grain == "c" & agreement.evAct$testCond == "sample"], agreement.evAct$sampleSize[agreement.evAct$grain == 'c' & agreement.evAct$testCond == "sample"], sampleSize)
agreement.fit.fine.evAct <- functionFit(agreement.evAct$agreementIndex[agreement.evAct$grain == "f" & agreement.evAct$testCond == "sample"], agreement.evAct$sampleSize[agreement.evAct$grain == 'f' & agreement.evAct$testCond == "sample"], sampleSize)
agreement.evAct.functionFit <- agreement.fit.coarse.evAct[[1]] %>% mutate(grain = "c")
agreement.evAct.functionFit <- rbind(agreement.evAct.functionFit, agreement.fit.fine.evAct[[1]] %>% mutate(grain = "f"))
agreement.evAct.functionFit$grain <- as.factor(agreement.evAct.functionFit$grain)
for(i in 1:nrow(agreement.evAct)){
if(agreement.evAct$grain[i] == 'c'){
agreement.evAct$elbow[i] = agreement.fit.coarse.evAct[[3]][1]
} else {
agreement.evAct$elbow[i] = agreement.fit.fine.evAct[[3]][1]
}
}Plot Agreement Index over sample size
Agreement Index over increasing sample sizes for segmentation of: A) commercial-lab and B) everyday-online.
agreement.comm.coarse.plot <- plotIndivAgreement(agreement.comm[agreement.comm$grain == 'c',], "agreementIndex", "Agreement Index", -0.0000000000000001, 1, "Coarse")
agreement.comm.fine.plot <- plotIndivAgreement(agreement.comm[agreement.comm$grain == 'f',],"agreementIndex", "Agreement Index", -0.0000000000000001, 1, "Fine")
agreement.evAct.coarse.plot <- plotIndivAgreement(agreement.evAct[agreement.evAct$grain == 'c', ], "agreementIndex", "Agreement Index", -0.0000000000000001, 1, "")
agreement.evAct.fine.plot <- plotIndivAgreement(agreement.evAct[agreement.evAct$grain == 'f', ], "agreementIndex", "Agreement Index", -0.0000000000000001, 1, "")
agreement.plot.combined <- ggarrange(agreement.comm.coarse.plot, agreement.comm.fine.plot,
agreement.evAct.coarse.plot, agreement.evAct.fine.plot,
nrow = 2, ncol = 2,
labels = c("A","", "B", ""),
common.legend = TRUE, legend = "bottom",
align = "hv")
agreement.plot.combined
setwd('../plots')
ggsave("esMethod_agreement.pdf", device = cairo_pdf, width = 25, height = 19, units = "cm", bg = "transparent")
# ggsave(plot = agreement.plot.combined, file = "esMethod_agreement.png",
# type = "cairo-png", bg = "transparent",
# width = 25, height = 19, units = "cm", dpi = 300)Plot Agreement Index function fits
Agreement Index fit to various growth and decay functions over increasing sample sizes for segmentation of: A) commercial-lab and B) everyday-online.
agreement.comm.coarse.funcPlot <- plotFunctionFits(agreement.functionFit[agreement.functionFit$grain == 'c',], agreement.comm[agreement.comm$grain == 'c' & agreement.comm$testCond == 'sample',], "agreementIndex", "Agreement Index", -0.0000000000000001, 1)
agreement.comm.fine.funcPlot <- plotFunctionFits(agreement.functionFit[agreement.functionFit$grain == 'f',], agreement.comm[agreement.comm$grain == 'f' & agreement.comm$testCond == 'sample',], "agreementIndex", "Agreement Index", -0.0000000000000001, 1)
agreement.evAct.coarse.funcPlot <- plotFunctionFits(agreement.evAct.functionFit[agreement.evAct.functionFit$grain == 'c',], agreement.evAct[agreement.evAct$grain == 'c' & agreement.evAct$testCond == 'sample',], "agreementIndex", "Agreement Index", -0.0000000000000001, 1)
agreement.evAct.fine.funcPlot <- plotFunctionFits(agreement.evAct.functionFit[agreement.evAct.functionFit$grain == 'f',], agreement.evAct[agreement.evAct$grain == 'f' & agreement.evAct$testCond == 'sample',], "agreementIndex", "Agreement Index", -0.0000000000000001, 1)
agreement.funcFit.plot.combined <- ggarrange(agreement.comm.coarse.funcPlot, agreement.comm.fine.funcPlot,
agreement.evAct.coarse.funcPlot, agreement.evAct.fine.funcPlot,
nrow = 2, ncol = 2,
labels = c("A","", "B", ""),
common.legend = TRUE, legend = "bottom",
align = "hv")
agreement.funcFit.plot.combined
setwd('../plots')
ggsave("esMethod_agreement_functionFit.pdf", device = cairo_pdf, width = 25, height = 19, units = "cm", bg = "transparent")Agreement Index linear models and analysis
#Commercial
agreement.comm.lmer <- lmer(agreementIndex~sampleSize*grain*testCond + (1|movieName), data = agreement.comm)
#emm_options(pbkrtest.limit = 3000)
agreement.comm.poly <- summary(contrast(emmeans(agreement.comm.lmer,'sampleSize', by = c( 'grain', 'testCond')), 'poly', max.degree = 3), infer = TRUE)#Polynomial contrast
agreement.comm.rand_vs_samp <- summary(emmeans(agreement.comm.lmer, c("testCond"), by = "grain", contr = list(c1 = c(0, -1, 1))), infer = T)# overall random - sample
agreement.comm.crossMov_vs_samp <- summary(emmeans(agreement.comm.lmer, c("testCond"), by = "grain", contr = list(c1 = c(-1, 0, 1))), infer = T)# overall crossMov - sample
agreement.comm.crossMov_vs_rand <- summary(emmeans(agreement.comm.lmer, c("testCond"), by = "grain", contr = list(c1 = c(-1, 1, 0))), infer = T)# overall crossMov - random
agreement.comm.pairwise <- summary(emmeans(agreement.comm.lmer, "testCond",by = c("sampleSize", "grain"), weights = "proportional", contr = "revpairwise", adjust = "none")$contrasts, infer = TRUE)#random vs sample pairwise/ sample size
agreement.evAct.lmer <- lmer(agreementIndex~sampleSize*grain*testCond + (1|movieName), data = agreement.evAct)
#emm_options(lmerTest.limit = 38688)
agreement.evAct.poly <- summary(contrast(emmeans(agreement.evAct.lmer,'sampleSize', by = c( 'grain', 'testCond')), 'poly', max.degree = 3), infer = TRUE)#Polynomial contrast
agreement.evAct.rand_vs_samp <- summary(emmeans(agreement.evAct.lmer, c("testCond"), by = "grain", contr = list(c1 = c(0, -1, 1))), infer = T)# overall random - sample
agreement.evAct.crossMov_vs_samp <- summary(emmeans(agreement.evAct.lmer, c("testCond"), by = "grain", contr = list(c1 = c(-1, 0, 1))), infer = T)# overall crossMov - sample
agreement.evAct.crossMov_vs_rand <- summary(emmeans(agreement.evAct.lmer, c("testCond"), by = "grain", contr = list(c1 = c(-1, 1, 0))), infer = T)# overall crossMov - random
agreement.evAct.pairwise <- summary(emmeans(agreement.evAct.lmer, "testCond",by = c("sampleSize", "grain"), weights = "proportional", contr = "revpairwise", adjust = "none")$contrasts, infer = TRUE)#random vs sample pairwise/ sample sizeAgreement Index random - sample contrast commercial-lab
agreement.comm.rand_vs_samp$contrasts## grain = c:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio p.value
## c1 0.1577242 0.001170359 Inf 0.1554303 0.1600180 134.766 <.0001
##
## grain = f:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio p.value
## c1 0.1205933 0.001170359 Inf 0.1182995 0.1228872 103.040 <.0001
##
## Results are averaged over the levels of: sampleSize
## Confidence level used: 0.95Agreement Index random - sample contrast everyday-online
agreement.evAct.rand_vs_samp$contrasts## grain = c:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio p.value
## c1 0.1904387 0.001051004 Inf 0.1883788 0.1924986 181.197 <.0001
##
## grain = f:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio p.value
## c1 0.1612430 0.001051004 Inf 0.1591831 0.1633029 153.418 <.0001
##
## Results are averaged over the levels of: sampleSize
## Confidence level used: 0.95Agreement Index cross movie - sample contrast commercial-lab
agreement.comm.crossMov_vs_samp$contrasts## grain = c:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio p.value
## c1 0.2291428 0.001170359 Inf 0.2268489 0.2314366 195.788 <.0001
##
## grain = f:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio p.value
## c1 0.2037672 0.001170359 Inf 0.2014733 0.2060611 174.107 <.0001
##
## Results are averaged over the levels of: sampleSize
## Confidence level used: 0.95Agreement Index cross movie - sample contrast everyday-online
agreement.evAct.crossMov_vs_samp$contrasts## grain = c:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio p.value
## c1 0.2232330 0.001051004 Inf 0.2211731 0.2252929 212.40 <.0001
##
## grain = f:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio p.value
## c1 0.1864058 0.001051004 Inf 0.1843458 0.1884657 177.36 <.0001
##
## Results are averaged over the levels of: sampleSize
## Confidence level used: 0.95Agreement Index cross movie - random contrast commercial-lab
agreement.comm.crossMov_vs_rand$contrasts## grain = c:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio p.value
## c1 0.07141858 0.001170359 Inf 0.06912472 0.07371244 61.023 <.0001
##
## grain = f:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio p.value
## c1 0.08317385 0.001170359 Inf 0.08087999 0.08546771 71.067 <.0001
##
## Results are averaged over the levels of: sampleSize
## Confidence level used: 0.95Agreement Index cross movie - random contrast everyday-online
agreement.evAct.crossMov_vs_rand$contrasts## grain = c:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio p.value
## c1 0.03279429 0.001051004 Inf 0.03073436 0.03485422 31.203 <.0001
##
## grain = f:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio p.value
## c1 0.02516278 0.001051004 Inf 0.02310285 0.02722271 23.942 <.0001
##
## Results are averaged over the levels of: sampleSize
## Confidence level used: 0.95Agreement Index smallest sample size for significant pairwise effect commercial-lab
agreement.comm.pairwise %>% dplyr::filter(p.value < .05, z.ratio > 0) %>% dplyr::group_by(grain, contrast) %>% filter(sampleSize == min(as.numeric(as.character(sampleSize))))## # A tibble: 6 x 10
## # Groups: grain, contrast [6]
## contrast sampleSize grain estimate SE df asymp.LCL asymp.UCL
## <fct> <fct> <fct> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 sample … 2 c 0.0566 0.00468 Inf 0.0474 0.0658
## 2 sample … 2 c 0.0898 0.00468 Inf 0.0806 0.0989
## 3 random … 6 c 0.0129 0.00468 Inf 0.00370 0.0220
## 4 sample … 2 f 0.0400 0.00468 Inf 0.0308 0.0491
## 5 sample … 2 f 0.101 0.00468 Inf 0.0918 0.110
## 6 random … 4 f 0.0117 0.00468 Inf 0.00248 0.0208
## # … with 2 more variables: z.ratio <dbl>, p.value <dbl>Agreement Index smallest sample size for significant pairwise effect everyday-online
agreement.evAct.pairwise %>% dplyr::filter(p.value < .05, z.ratio > 0) %>% dplyr::group_by(grain, contrast) %>% filter(sampleSize == min(as.numeric(as.character(sampleSize))))## # A tibble: 6 x 10
## # Groups: grain, contrast [6]
## contrast sampleSize grain estimate SE df asymp.LCL asymp.UCL
## <fct> <fct> <fct> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 sample … 2 c 0.0469 0.00420 Inf 0.0386 0.0551
## 2 sample … 2 c 0.0803 0.00420 Inf 0.0721 0.0886
## 3 random … 12 c 0.0169 0.00420 Inf 0.00862 0.0251
## 4 sample … 2 f 0.0659 0.00420 Inf 0.0576 0.0741
## 5 sample … 2 f 0.134 0.00420 Inf 0.126 0.142
## 6 random … 6 f 0.0109 0.00420 Inf 0.00262 0.0191
## # … with 2 more variables: z.ratio <dbl>, p.value <dbl>Agreement Index smallest pairwise effect commercial-lab
agreement.comm.pairwise %>% dplyr::filter(p.value < .05, z.ratio > 0) %>% dplyr::group_by(grain, contrast) %>% filter(z.ratio == min(z.ratio))## # A tibble: 6 x 10
## # Groups: grain, contrast [6]
## contrast sampleSize grain estimate SE df asymp.LCL asymp.UCL
## <fct> <fct> <fct> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 sample … 2 c 0.0566 0.00468 Inf 0.0474 0.0658
## 2 sample … 2 c 0.0898 0.00468 Inf 0.0806 0.0989
## 3 random … 6 c 0.0129 0.00468 Inf 0.00370 0.0220
## 4 sample … 2 f 0.0400 0.00468 Inf 0.0308 0.0491
## 5 sample … 2 f 0.101 0.00468 Inf 0.0918 0.110
## 6 random … 4 f 0.0117 0.00468 Inf 0.00248 0.0208
## # … with 2 more variables: z.ratio <dbl>, p.value <dbl>Agreement Index smallest pairwise effect everyday-online
agreement.evAct.pairwise %>% dplyr::filter(p.value < .05, z.ratio > 0) %>% dplyr::group_by(grain, contrast) %>% filter(z.ratio == min(z.ratio))## # A tibble: 6 x 10
## # Groups: grain, contrast [6]
## contrast sampleSize grain estimate SE df asymp.LCL asymp.UCL
## <fct> <fct> <fct> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 sample … 2 c 0.0469 0.00420 Inf 0.0386 0.0551
## 2 sample … 2 c 0.0803 0.00420 Inf 0.0721 0.0886
## 3 random … 12 c 0.0169 0.00420 Inf 0.00862 0.0251
## 4 sample … 2 f 0.0659 0.00420 Inf 0.0576 0.0741
## 5 sample … 2 f 0.134 0.00420 Inf 0.126 0.142
## 6 random … 6 f 0.0109 0.00420 Inf 0.00262 0.0191
## # … with 2 more variables: z.ratio <dbl>, p.value <dbl>Agreement Index polynomial fit commercial-lab
agreement.comm.poly[agreement.comm.poly$testCond == "sample",]## contrast grain testCond estimate SE df asymp.LCL asymp.UCL
## 13 linear c sample 10.18298 0.1220769 Inf 9.94371 10.422243
## 14 quadratic c sample -10.15537 0.2501832 Inf -10.64572 -9.665016
## 15 cubic c sample 71.02177 3.3230932 Inf 64.50863 77.534917
## 16 linear f sample 10.38112 0.1220769 Inf 10.14185 10.620384
## 17 quadratic f sample -11.83477 0.2501832 Inf -12.32512 -11.344421
## 18 cubic f sample 90.35719 3.3230932 Inf 83.84405 96.870336
## z.ratio p.value
## 13 83.41446 0.000000e+00
## 14 -40.59173 0.000000e+00
## 15 21.37219 2.424961e-101
## 16 85.03754 0.000000e+00
## 17 -47.30443 0.000000e+00
## 18 27.19069 8.369360e-163Agreement Index polynomial fit everyday-online
agreement.evAct.poly[agreement.evAct.poly$testCond == "sample",]## contrast grain testCond estimate SE df asymp.LCL asymp.UCL
## 13 linear c sample 11.456484 0.1096273 Inf 11.241618 11.671350
## 14 quadratic c sample -10.784017 0.2246692 Inf -11.224360 -10.343673
## 15 cubic c sample 70.339801 2.9842000 Inf 64.490877 76.188726
## 16 linear f sample 8.101342 0.1096273 Inf 7.886476 8.316208
## 17 quadratic f sample -9.983328 0.2246692 Inf -10.423671 -9.542984
## 18 cubic f sample 78.066083 2.9842000 Inf 72.217158 83.915007
## z.ratio p.value
## 13 104.50391 0.000000e+00
## 14 -47.99954 0.000000e+00
## 15 23.57074 7.694444e-123
## 16 73.89893 0.000000e+00
## 17 -44.43568 0.000000e+00
## 18 26.15980 7.623629e-151Surprise Index
Load dataset
setwd('../data/bootstrapped')
surprise.comm <- read.delim("esMethods_SurpriseIndex_c0.1_f0.05_100Iterations_commercial.txt", head = TRUE)
surprise.comm$sampleSize <- as.factor(surprise.comm$sampleSize)
surprise.comm$logSurprise <- log(surprise.comm$surpriseIndex, base = 10)
surprise.evAct <- read.delim("esMethods_SurpriseIndex_c0.1_f0.05_100Iterations_evAct.txt", head = TRUE)
surprise.evAct$sampleSize <- as.factor(surprise.evAct$sampleSize)
surprise.evAct$logSurprise <- log(surprise.evAct$surpriseIndex, base = 10)Fit Surprise Index values to growth and decay functions
##commercial
surprise.fit.coarse.comm <- functionFit(surprise.comm$surpriseIndex[surprise.comm$grain == "c" & surprise.comm$testCond == 'sample'], surprise.comm$sampleSize[surprise.comm$grain == "c" & surprise.comm$testCond == 'sample'], sampleSize)
surprise.fit.fine.comm <- functionFit(surprise.comm$surpriseIndex[surprise.comm$grain == "f" & surprise.comm$testCond == 'sample'], surprise.comm$sampleSize[surprise.comm$grain == "f" & surprise.comm$testCond == 'sample'], sampleSize)
surprise.functionFit <- surprise.fit.coarse.comm[[1]] %>% mutate(grain = "c")
surprise.functionFit <- rbind(surprise.functionFit, surprise.fit.fine.comm[[1]] %>% mutate(grain = "f"))
surprise.functionFit$grain <- as.factor(surprise.functionFit$grain)
#surprise.functionFit.comm.figure <- plotFunctionFits(surprise.functionFit, surprise.comm[surprise.comm$testCond == "sample",], "surpriseIndex", "Surprise Index (raw)")
for(i in 1:nrow(surprise.comm)){
if(surprise.comm$grain[i] == 'c'){
surprise.comm$elbow[i] = surprise.fit.coarse.comm[[3]][1]
} else {
surprise.comm$elbow[i] = surprise.fit.fine.comm[[3]][1]
}
}
##everyday Activities
surprise.fit.coarse.evAct <- functionFit(surprise.evAct$surpriseIndex[surprise.evAct$grain == "c" & surprise.evAct$testCond == "sample"], surprise.evAct$sampleSize[surprise.evAct$grain == 'c' & surprise.evAct$testCond == "sample"], sampleSize)
surprise.fit.fine.evAct <- functionFit(surprise.evAct$surpriseIndex[surprise.evAct$grain == "f" & surprise.evAct$testCond == "sample"], surprise.evAct$sampleSize[surprise.evAct$grain == 'f' & surprise.evAct$testCond == "sample"], sampleSize)
surprise.evAct.functionFit <- surprise.fit.coarse.evAct[[1]] %>% mutate(grain = "c")
surprise.evAct.functionFit <- rbind(surprise.evAct.functionFit, surprise.fit.fine.evAct[[1]] %>% mutate(grain = "f"))
surprise.evAct.functionFit$grain <- as.factor(surprise.evAct.functionFit$grain)
#surprise.functionFit.evAct.figure <- plotFunctionFits(surprise.evAct.functionFit, surprise.evAct[surprise.evAct$testCond == "sample",], "surpriseIndex", "Surprise Index (raw)") + ylim(0,15)
for(i in 1:nrow(surprise.evAct)){
if(surprise.evAct$grain[i] == 'c'){
surprise.evAct$elbow[i] = surprise.fit.coarse.evAct[[3]][1]
} else {
surprise.evAct$elbow[i] = surprise.fit.fine.evAct[[3]][1]
}
}Plot log10 Surprise Index over sample size
Log10 transformed peak-to-peak distance over increasing sample sizes for segmentation of: A) commercial-lab and B) everyday-online.
surprise.comm.coarse.plot <- plotIndivAgreement(surprise.comm[surprise.comm$grain == 'c',], "logSurprise", "Surprise Index (log)", -1, 2, "Coarse")
surprise.comm.fine.plot <- plotIndivAgreement(surprise.comm[surprise.comm$grain == 'f',], "logSurprise", "Surprise Index (log)", -1, 2, "Fine")
surprise.evAct.coarse.plot <- plotIndivAgreement(surprise.evAct[surprise.evAct$grain == 'c', ], "logSurprise", "Surprise Index (log)", -1, 2, "")
surprise.evAct.fine.plot <- plotIndivAgreement(surprise.evAct[surprise.evAct$grain == 'f', ], "logSurprise", "Surprise Index (log)", -1, 2, "")
surprise.plot.combined <- ggarrange(surprise.comm.coarse.plot, surprise.comm.fine.plot,
surprise.evAct.coarse.plot, surprise.evAct.fine.plot,
nrow = 2, ncol = 2,
labels = c("A","", "B", ""),
common.legend = TRUE, legend = "bottom",
align = "hv")
surprise.plot.combined
setwd('../plots')
ggsave("esMethod_surprise.pdf", device = cairo_pdf, width = 25, height = 19, units = "cm", bg = "transparent", dpi = 300)Plot Surprise Index function fits
Peak-to-peak values fit to various growth and decay functions over increasing sample sizes for segmentation of: A) commercial-lab and B) everyday-online.
surprise.comm.coarse.funcPlot <- plotFunctionFits(surprise.functionFit[surprise.functionFit$grain == 'c',], surprise.comm[surprise.comm$grain == 'c' & surprise.comm$testCond == 'sample',], "surpriseIndex", "Surprise Index", 0, 15)
surprise.comm.fine.funcPlot <- plotFunctionFits(surprise.functionFit[surprise.functionFit$grain == 'f',], surprise.comm[surprise.comm$grain == 'f' & surprise.comm$testCond == 'sample',], "surpriseIndex", "Surprise Index", 0, 15)
surprise.evAct.coarse.funcPlot <- plotFunctionFits(surprise.evAct.functionFit[surprise.evAct.functionFit$grain == 'c',], surprise.evAct[surprise.evAct$grain == 'c' & surprise.evAct$testCond == 'sample',], "surpriseIndex", "Surprise Index", 0, 15)
surprise.evAct.fine.funcPlot <- plotFunctionFits(surprise.evAct.functionFit[surprise.evAct.functionFit$grain == 'f',], surprise.evAct[surprise.evAct$grain == 'f' & surprise.evAct$testCond == 'sample',], "surpriseIndex", "Surprise Index", 0, 15)
surprise.funcFit.plot.combined <- ggarrange(surprise.comm.coarse.funcPlot, surprise.comm.fine.funcPlot,
surprise.evAct.coarse.funcPlot, surprise.evAct.fine.funcPlot,
nrow = 2, ncol = 2,
labels = c("A","", "B", ""),
common.legend = TRUE, legend = "bottom",
align = "hv")
surprise.funcFit.plot.combined
setwd('../plots')
ggsave("esMethod_surprise_functionFit.pdf", device = cairo_pdf, width = 25, height = 19, units = "cm", bg = "transparent")Surprise Index linear models and analysis
#Commercial
surprise.comm.lmer <- lmer(surpriseIndex~sampleSize*grain*testCond + (1|movieName), data = surprise.comm)
#emm_options(pbkrtest.limit = 3000)
surprise.comm.poly <- summary(contrast(emmeans(surprise.comm.lmer,'sampleSize', by = c( 'grain', 'testCond')), 'poly', max.degree = 3), infer = TRUE)#Polynomial contrast
surprise.comm.rand_vs_samp <- summary(emmeans(surprise.comm.lmer, c("testCond"), by = "grain", contr = list(c1 = c(0, -1, 1))), infer = T)# overall random - sample
surprise.comm.crossMov_vs_samp <- summary(emmeans(surprise.comm.lmer, c("testCond"), by = "grain", contr = list(c1 = c(-1, 0, 1))), infer = T)# overall crossMov - sample
surprise.comm.crossMov_vs_rand <- summary(emmeans(surprise.comm.lmer, c("testCond"), by = "grain", contr = list(c1 = c(-1, 1, 0))), infer = T)# overall crossMov - random
surprise.comm.pairwise <- summary(emmeans(surprise.comm.lmer, "testCond",by = c("sampleSize", "grain"), weights = "proportional", contr = "revpairwise", adjust = "none")$contrasts, infer = TRUE)#random vs sample pairwise/ sample size
surprise.evAct.lmer <- lmer(surpriseIndex~sampleSize*grain*testCond + (1|movieName), data = surprise.evAct)
#emm_options(lmerTest.limit = 38400)
surprise.evAct.poly <- summary(contrast(emmeans(surprise.evAct.lmer,'sampleSize', by = c( 'grain', 'testCond')), 'poly', max.degree = 3), infer = TRUE)#Polynomial contrast
surprise.evAct.rand_vs_samp <- summary(emmeans(surprise.evAct.lmer, c("testCond"), by = "grain", contr = list(c1 = c(0, -1, 1))), infer = T)# overall random - sample
surprise.evAct.crossMov_vs_samp <- summary(emmeans(surprise.evAct.lmer, c("testCond"), by = "grain", contr = list(c1 = c(-1, 0, 1))), infer = T)# overall crossMov - sample
surprise.evAct.crossMov_vs_rand <- summary(emmeans(surprise.evAct.lmer, c("testCond"), by = "grain", contr = list(c1 = c(1, -1, 0))), infer = T)# overall crossMov - random
surprise.evAct.pairwise <- summary(emmeans(surprise.evAct.lmer, "testCond",by = c("sampleSize", "grain"), weights = "proportional", contr = "revpairwise", adjust = "none")$contrasts, infer = TRUE)#random vs sample pairwise/ sample sizeSurprise Index random - sample contrast commercial-lab
surprise.comm.rand_vs_samp$contrasts## grain = c:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio p.value
## c1 3.434135 0.03548559 Inf 3.364585 3.503686 96.775 <.0001
##
## grain = f:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio p.value
## c1 6.705420 0.03548559 Inf 6.635870 6.774971 188.962 <.0001
##
## Results are averaged over the levels of: sampleSize
## Confidence level used: 0.95Surprise Index random - sample contrast everyday-online
surprise.evAct.rand_vs_samp$contrasts## grain = c:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio p.value
## c1 2.549908 0.02249992 Inf 2.505809 2.594007 113.330 <.0001
##
## grain = f:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio p.value
## c1 4.052439 0.02249992 Inf 4.008340 4.096538 180.109 <.0001
##
## Results are averaged over the levels of: sampleSize
## Confidence level used: 0.95Surprise Index cross movie - sample contrast commercial-lab
surprise.comm.crossMov_vs_samp$contrasts## grain = c:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio p.value
## c1 3.570111 0.03548559 Inf 3.500561 3.639662 100.607 <.0001
##
## grain = f:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio p.value
## c1 6.679965 0.03548559 Inf 6.610415 6.749516 188.244 <.0001
##
## Results are averaged over the levels of: sampleSize
## Confidence level used: 0.95Surprise Index cross movie - sample contrast everyday-online
surprise.evAct.crossMov_vs_samp$contrasts## grain = c:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio p.value
## c1 2.262476 0.02249992 Inf 2.218377 2.306575 100.555 <.0001
##
## grain = f:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio p.value
## c1 3.885204 0.02249992 Inf 3.841105 3.929303 172.676 <.0001
##
## Results are averaged over the levels of: sampleSize
## Confidence level used: 0.95Surprise Index cross movie - random contrast commercial-lab
surprise.comm.crossMov_vs_rand$contrasts## grain = c:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio
## c1 0.13597635 0.03548559 Inf 0.06642588 0.20552683 3.832
## p.value
## 0.0001
##
## grain = f:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio
## c1 -0.02545493 0.03548559 Inf -0.09500541 0.04409554 -0.717
## p.value
## 0.4732
##
## Results are averaged over the levels of: sampleSize
## Confidence level used: 0.95Surprise Index cross movie - random contrast everyday-online
surprise.evAct.crossMov_vs_rand$contrasts## grain = c:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio p.value
## c1 0.2874323 0.02249992 Inf 0.2433333 0.3315313 12.775 <.0001
##
## grain = f:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio p.value
## c1 0.1672345 0.02249992 Inf 0.1231355 0.2113335 7.433 <.0001
##
## Results are averaged over the levels of: sampleSize
## Confidence level used: 0.95Surprise Index smallest sample size for significant pairwise effect commercial-lab
surprise.comm.pairwise %>% dplyr::filter(p.value < .05, z.ratio > 0) %>% dplyr::group_by(grain, contrast) %>% filter(sampleSize == min(as.numeric(as.character(sampleSize))))## # A tibble: 4 x 10
## # Groups: grain, contrast [4]
## contrast sampleSize grain estimate SE df asymp.LCL asymp.UCL
## <fct> <fct> <fct> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 sample … 2 c 1.77 0.142 Inf 1.49 2.04
## 2 sample … 2 c 1.60 0.142 Inf 1.33 1.88
## 3 sample … 2 f 2.33 0.142 Inf 2.05 2.60
## 4 sample … 2 f 2.32 0.142 Inf 2.04 2.60
## # … with 2 more variables: z.ratio <dbl>, p.value <dbl>Surprise Index smallest sample size for significant pairwise effect everyday-online
surprise.evAct.pairwise %>% dplyr::filter(p.value < .05, z.ratio > 0) %>% dplyr::group_by(grain, contrast) %>% filter(sampleSize == min(as.numeric(as.character(sampleSize))))## # A tibble: 4 x 10
## # Groups: grain, contrast [4]
## contrast sampleSize grain estimate SE df asymp.LCL asymp.UCL
## <fct> <fct> <fct> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 sample … 2 c 1.28 0.0900 Inf 1.11 1.46
## 2 sample … 2 c 1.21 0.0900 Inf 1.03 1.38
## 3 sample … 2 f 1.88 0.0900 Inf 1.70 2.05
## 4 sample … 2 f 1.91 0.0900 Inf 1.73 2.09
## # … with 2 more variables: z.ratio <dbl>, p.value <dbl>Surprise Index smallest pairwise effect commercial-lab
surprise.comm.pairwise %>% dplyr::filter(p.value < .05, z.ratio > 0) %>% dplyr::group_by(grain, contrast) %>% filter(z.ratio == min(z.ratio))## # A tibble: 4 x 10
## # Groups: grain, contrast [4]
## contrast sampleSize grain estimate SE df asymp.LCL asymp.UCL
## <fct> <fct> <fct> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 sample … 4 c 1.44 0.142 Inf 1.16 1.72
## 2 sample … 4 c 1.33 0.142 Inf 1.05 1.61
## 3 sample … 2 f 2.33 0.142 Inf 2.05 2.60
## 4 sample … 2 f 2.32 0.142 Inf 2.04 2.60
## # … with 2 more variables: z.ratio <dbl>, p.value <dbl>Surprise Index smallest pairwise effect everyday-online
surprise.evAct.pairwise %>% dplyr::filter(p.value < .05, z.ratio > 0) %>% dplyr::group_by(grain, contrast) %>% filter(z.ratio == min(z.ratio))## # A tibble: 4 x 10
## # Groups: grain, contrast [4]
## contrast sampleSize grain estimate SE df asymp.LCL asymp.UCL
## <fct> <fct> <fct> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 sample … 2 c 1.28 0.0900 Inf 1.11 1.46
## 2 sample … 2 c 1.21 0.0900 Inf 1.03 1.38
## 3 sample … 2 f 1.88 0.0900 Inf 1.70 2.05
## 4 sample … 2 f 1.91 0.0900 Inf 1.73 2.09
## # … with 2 more variables: z.ratio <dbl>, p.value <dbl>Surprise Index polynomial fit commercial-lab
surprise.comm.poly[surprise.comm.poly$testCond == "sample",]## contrast grain testCond estimate SE df asymp.LCL asymp.UCL
## 13 linear c sample 157.21190 3.701404 Inf 149.9573 164.46652
## 14 quadratic c sample -43.69346 7.585620 Inf -58.5610 -28.82591
## 15 cubic c sample -392.21574 100.757076 Inf -589.6960 -194.73550
## 16 linear f sample 287.41961 3.701404 Inf 280.1650 294.67423
## 17 quadratic f sample -211.13912 7.585620 Inf -226.0067 -196.27157
## 18 cubic f sample 351.69993 100.757076 Inf 154.2197 549.18017
## z.ratio p.value
## 13 42.473591 0.000000e+00
## 14 -5.760038 8.409524e-09
## 15 -3.892687 9.914009e-05
## 16 77.651517 0.000000e+00
## 17 -27.834127 1.676671e-170
## 18 3.490573 4.819859e-04Surprise Index polynomial fit everyday-online
surprise.evAct.poly[surprise.evAct.poly$testCond == "sample",]## contrast grain testCond estimate SE df asymp.LCL asymp.UCL
## 13 linear c sample 78.15937 2.346904 Inf 73.55952 82.75922
## 14 quadratic c sample -31.27365 4.809722 Inf -40.70053 -21.84676
## 15 cubic c sample -74.87400 63.885815 Inf -200.08790 50.33990
## 16 linear f sample 97.38823 2.346904 Inf 92.78838 101.98807
## 17 quadratic f sample -85.25332 4.809722 Inf -94.68021 -75.82644
## 18 cubic f sample 577.77115 63.885815 Inf 452.55725 702.98505
## z.ratio p.value
## 13 33.303180 3.472692e-243
## 14 -6.502173 7.916788e-11
## 15 -1.171997 2.411982e-01
## 16 41.496467 0.000000e+00
## 17 -17.725208 2.679205e-70
## 18 9.043810 1.513047e-19Predictive value and detection accuracy
Load dataset
setwd('../data/bootstrapped')
normativeHit.comm <- read.delim("esMethods_normativeHitRate_bin_100Iterations_commercial.txt", head = TRUE)
predictive.value<- normativeHit.comm %>% dplyr::select(Iteration, sampleSize, movieName, grain, testCond, normativePPV, normativeFalsePPV)
predictive.value$dPrime <- qnorm(predictive.value$normativePPV) - qnorm(predictive.value$normativeFalsePPV)
predictive.value <- predictive.value %>% dplyr::filter(is.finite(predictive.value$dPrime))
predictive.value$sampleSize <- as.factor(predictive.value$sampleSize)
detection.sensitivity <- normativeHit.comm %>% dplyr::select(Iteration, sampleSize, movieName, grain, testCond, normativeHit, normativeFA, normativeMiss, normativeCR)
detection.sensitivity$dPrime <- qnorm(detection.sensitivity$normativeHit) - qnorm(detection.sensitivity$normativeFA)
detection.sensitivity <- detection.sensitivity %>% dplyr::filter(is.finite(detection.sensitivity$dPrime))
detection.sensitivity$sampleSize <- as.factor(detection.sensitivity$sampleSize)Fit predictive value and detection accuracy values to growth and decay functions
#zalladPrime
predictive.fit.coarse <- functionFit(predictive.value$dPrime[predictive.value$grain == "c" & predictive.value$testCond == 'sample'], predictive.value$sampleSize[predictive.value$grain == "c" & predictive.value$testCond == 'sample'], sampleSize)## Error in optim(startVec, opfct, hessian = TRUE, method = optMethod, control = list(maxit = maxIt, :
## non-finite value supplied by optim
## Error in optim(startVec, opfct, hessian = TRUE, method = optMethod, control = list(maxit = maxIt, :
## non-finite value supplied by optimpredictive.fit.coarse <- functionFit(predictive.value$dPrime[predictive.value$grain == "f" & predictive.value$testCond == 'sample'], predictive.value$sampleSize[predictive.value$grain == "f" & predictive.value$testCond == 'sample'], sampleSize)## Error in optim(startVec, opfct, hessian = TRUE, method = optMethod, control = list(maxit = maxIt, :
## non-finite value supplied by optim
## Error in optim(startVec, opfct, hessian = TRUE, method = optMethod, control = list(maxit = maxIt, :
## non-finite value supplied by optimpredictive.functionFit <- predictive.fit.coarse[[1]] %>% mutate(grain = "c")
predictive.functionFit <- rbind(predictive.functionFit, predictive.fit.coarse[[1]] %>% mutate(grain = "f"))
predictive.functionFit$grain <- as.factor(predictive.functionFit$grain)
for(i in 1:nrow(predictive.value)){
if(predictive.value$grain[i] == 'c'){
predictive.value$elbow[i] = predictive.fit.coarse[[3]][1]
} else {
predictive.value$elbow[i] = predictive.fit.coarse[[3]][1]
}
}
##commercial
detection.sensitivity.fit.coarse.comm <- functionFit(detection.sensitivity$dPrime[detection.sensitivity$grain == "c" & detection.sensitivity$testCond == 'sample'], detection.sensitivity$sampleSize[detection.sensitivity$grain == "c" & detection.sensitivity$testCond == 'sample'], sampleSize)
detection.sensitivity.fit.fine.comm <- functionFit(detection.sensitivity$dPrime[detection.sensitivity$grain == "f" & detection.sensitivity$testCond == 'sample'], detection.sensitivity$sampleSize[detection.sensitivity$grain == "f" & detection.sensitivity$testCond == 'sample'], sampleSize)
detection.sensitivity.functionFit <- detection.sensitivity.fit.coarse.comm[[1]] %>% mutate(grain = "c")
detection.sensitivity.functionFit <- rbind(detection.sensitivity.functionFit, detection.sensitivity.fit.fine.comm[[1]] %>% mutate(grain = "f"))
detection.sensitivity.functionFit$grain <- as.factor(detection.sensitivity.functionFit$grain)
for(i in 1:nrow(detection.sensitivity)){
if(detection.sensitivity$grain[i] == 'c'){
detection.sensitivity$elbow[i] = detection.sensitivity.fit.coarse.comm[[3]][1]
} else {
detection.sensitivity$elbow[i] = detection.sensitivity.fit.fine.comm[[3]][1]
}
}Plot predictive value and detection accuracy over sample size
Predictive value and detection accuracy over increasing sample sizes for segmentation of: A) commercial-lab and B) everyday-online.
predictive.value.comm.coarse.plot <- plotIndivAgreement(predictive.value[predictive.value$grain == 'c',], "dPrime", "Predictive value", -6, 1, "Coarse")
predictive.value.comm.fine.plot <- plotIndivAgreement(predictive.value[predictive.value$grain == 'f',],"dPrime", "Predictive value", -6, 1, "Fine")
predictive.value.plot <- ggarrange(predictive.value.comm.coarse.plot, predictive.value.comm.fine.plot,
nrow = 1, ncol = 2,
common.legend = TRUE, legend = "bottom",
align = "hv")
detection.sensitivity.coarse.plot <- plotIndivAgreement(detection.sensitivity[detection.sensitivity$grain == 'c',], "dPrime", "Detection accuracy", -1, 2, "Coarse")
detection.sensitivity.fine.plot <- plotIndivAgreement(detection.sensitivity[detection.sensitivity$grain == 'f',], "dPrime", "Detection accuracy", -1, 2, "Fine")
detection.sensitivity.plot <- ggarrange(detection.sensitivity.coarse.plot, detection.sensitivity.fine.plot,
nrow = 1, ncol = 2,
common.legend = TRUE, legend = "bottom",
align = "hv")
predictive.detection.accuracy <- ggarrange(predictive.value.plot, detection.sensitivity.plot,
nrow = 2, ncol = 1,
common.legend = TRUE, legend = "bottom",
labels = c("A","B"),
align = "hv")
predictive.detection.accuracy
setwd('../plots')
ggsave("esMethod_predictive&detection.pdf", device = cairo_pdf, width = 25, height = 19, units = "cm", bg = "transparent")Plot predictive value and detection accuracy function fits
Predictive value and detection accuracy fit to various growth and decay functions over increasing sample sizes for segmentation of: A) commercial-lab and B) everyday-online.
predictive.coarse.funcPlot <- plotFunctionFits(predictive.functionFit[predictive.functionFit$grain == 'c',], predictive.value[predictive.value$grain == 'c' & predictive.value$testCond == 'sample',], "dPrime", "Predictive value", -4, 0)
predictive.fine.funcPlot <- plotFunctionFits(predictive.functionFit[predictive.functionFit$grain == 'f',], predictive.value[predictive.value$grain == 'f' & predictive.value$testCond == 'sample',], "dPrime", "Predictive value", -4, 0)
detection.coarse.funcPlot <- plotFunctionFits(detection.sensitivity.functionFit[detection.sensitivity.functionFit$grain == 'c',], detection.sensitivity[detection.sensitivity$grain == 'c' & detection.sensitivity$testCond == 'sample',], "dPrime", "Detection accuracy", -1, 2)
detection.fine.funcPlot <- plotFunctionFits(detection.sensitivity.functionFit[detection.sensitivity.functionFit$grain == 'f',], detection.sensitivity[detection.sensitivity$grain == 'f' & detection.sensitivity$testCond == 'sample',], "dPrime", "Detection accuracy", -1, 2)
predictive.detection.funcFit.plot.combined <- ggarrange(predictive.coarse.funcPlot, predictive.fine.funcPlot,
detection.coarse.funcPlot, detection.fine.funcPlot,
nrow = 2, ncol = 2,
labels = c("A","", "B", ""),
common.legend = TRUE, legend = "bottom",
align = "hv")
predictive.detection.funcFit.plot.combined
setwd('../plots')
ggsave("esMethod_predictive&detection_functionFit.pdf", device = cairo_pdf, width = 25, height = 19, units = "cm", bg = "transparent")Predictive value and detection accuracy linear models and analysis
#Predictive value
predictive.value.lmer <- lmer(dPrime~sampleSize*grain*testCond + (1|movieName), data = predictive.value)
#emm_options(pbkrtest.limit = 3000)
predictive.rand_vs_samp <- summary(emmeans(predictive.value.lmer, c("testCond"), by = "grain", contr = list(c1 = c(0, -1, 1))), infer = T)# overall random - sample
predictive.crossMov_vs_samp <- summary(emmeans(predictive.value.lmer, c("testCond"), by = "grain", contr = list(c1 = c(-1, 0, 1))), infer = T)# overall crossMov - sample
predictive.crossMov_vs_rand <- summary(emmeans(predictive.value.lmer, c("testCond"), by = "grain", contr = list(c1 = c(1, -1, 0))), infer = T)# overall crossMov - random
predictive.poly <- summary(contrast(emmeans(predictive.value.lmer,'sampleSize', by = c( 'grain', 'testCond')), 'poly', max.degree = 3), infer = TRUE)#Polynomial contrast
predictive.pairwise <- summary(emmeans(predictive.value.lmer, "testCond",by = c("sampleSize", "grain"), weights = "proportional", contr = "revpairwise", adjust = "none")$contrasts, infer = TRUE)#random vs sample pairwise/ sample size
#Detection sensitivity
detection.sensitivity.lmer <- lmer(dPrime~sampleSize*grain*testCond + (1|movieName), data = detection.sensitivity)
#emm_options(pbkrtest.limit = 19029)
detection.poly <- summary(contrast(emmeans(detection.sensitivity.lmer,'sampleSize', by = c( 'grain', 'testCond')), 'poly', max.degree = 3), infer = TRUE)#Polynomial contrast
detection.rand_vs_samp <- summary(emmeans(detection.sensitivity.lmer, c("testCond"), by = "grain", contr = list(c1 = c(0, -1, 1))), infer = T)# overall random - sample
detection.crossMov_vs_samp <- summary(emmeans(detection.sensitivity.lmer, c("testCond"), by = "grain", contr = list(c1 = c(-1, 0, 1))), infer = T)# overall crossMov - sample
detection.crossMov_vs_rand <- summary(emmeans(detection.sensitivity.lmer, c("testCond"), by = "grain", contr = list(c1 = c(1, -1, 0))), infer = T)# overall crossMov - random
detection.pairwise <- summary(emmeans(detection.sensitivity.lmer, "testCond",by = c("sampleSize", "grain"), weights = "proportional", contr = "revpairwise", adjust = "none")$contrasts, infer = TRUE)#random vs sample pairwise/ sample sizePredictive value random - sample contrast commercial-lab
predictive.rand_vs_samp$contrasts## grain = c:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio p.value
## c1 1.585063 0.009201382 Inf 1.567029 1.603097 172.264 <.0001
##
## grain = f:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio p.value
## c1 1.115986 0.009081682 Inf 1.098186 1.133786 122.883 <.0001
##
## Results are averaged over the levels of: sampleSize
## Confidence level used: 0.95Detection value random - sample contrast commercial-lab
detection.rand_vs_samp$contrasts## grain = c:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio p.value
## c1 0.8096764 0.004325468 Inf 0.8011986 0.8181541 187.188 <.0001
##
## grain = f:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio p.value
## c1 0.6529486 0.004269173 Inf 0.6445812 0.6613161 152.945 <.0001
##
## Results are averaged over the levels of: sampleSize
## Confidence level used: 0.95Predictive value cross movie - sample contrast commercial-lab
predictive.crossMov_vs_samp$contrasts## grain = c:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio p.value
## c1 1.481705 0.009156284 Inf 1.463759 1.499651 161.824 <.0001
##
## grain = f:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio p.value
## c1 1.125738 0.009081682 Inf 1.107938 1.143538 123.957 <.0001
##
## Results are averaged over the levels of: sampleSize
## Confidence level used: 0.95Detection value cross movie - sample contrast commercial-lab
detection.crossMov_vs_samp$contrasts## grain = c:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio p.value
## c1 0.7738612 0.004304249 Inf 0.7654250 0.7822974 179.790 <.0001
##
## grain = f:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio p.value
## c1 0.6622736 0.004269173 Inf 0.6539062 0.6706410 155.129 <.0001
##
## Results are averaged over the levels of: sampleSize
## Confidence level used: 0.95Predictive value cross movie - random contrast commercial-lab
predictive.crossMov_vs_rand$contrasts## grain = c:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio
## c1 0.103357684 0.009176672 Inf 0.08537174 0.121343630 11.263
## p.value
## <.0001
##
## grain = f:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio
## c1 -0.009752199 0.009079522 Inf -0.02754774 0.008043337 -1.074
## p.value
## 0.2828
##
## Results are averaged over the levels of: sampleSize
## Confidence level used: 0.95Detection value cross movie - random contrast commercial-lab
detection.crossMov_vs_rand$contrasts## grain = c:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio
## c1 0.03581518 0.004313832 Inf 0.02736022 0.0442701338 8.302
## p.value
## <.0001
##
## grain = f:
## contrast estimate SE df asymp.LCL asymp.UCL z.ratio
## c1 -0.00932495 0.004268157 Inf -0.01769038 -0.0009595155 -2.185
## p.value
## 0.0289
##
## Results are averaged over the levels of: sampleSize
## Confidence level used: 0.95Predictive value smallest sample size for significant pairwise effect commercial-lab
predictive.pairwise %>% dplyr::filter(p.value < .05, z.ratio > 0) %>% dplyr::group_by(grain, contrast) %>% filter(sampleSize == min(as.numeric(as.character(sampleSize))))## # A tibble: 5 x 10
## # Groups: grain, contrast [5]
## contrast sampleSize grain estimate SE df asymp.LCL asymp.UCL
## <fct> <fct> <fct> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 sample … 2 c 0.737 0.0398 Inf 0.659 0.815
## 2 sample … 2 c 0.914 0.0422 Inf 0.831 0.997
## 3 random … 4 c 0.0978 0.0377 Inf 0.0240 0.172
## 4 sample … 2 f 0.439 0.0365 Inf 0.368 0.511
## 5 sample … 2 f 0.516 0.0365 Inf 0.444 0.587
## # … with 2 more variables: z.ratio <dbl>, p.value <dbl>Detection value smallest sample size for significant pairwise effect commercial-lab
detection.pairwise %>% dplyr::filter(p.value < .05, z.ratio > 0) %>% dplyr::group_by(grain, contrast) %>% filter(sampleSize == min(as.numeric(as.character(sampleSize))))## # A tibble: 6 x 10
## # Groups: grain, contrast [6]
## contrast sampleSize grain estimate SE df asymp.LCL asymp.UCL
## <fct> <fct> <fct> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 sample … 2 c 0.343 0.0187 Inf 0.306 0.380
## 2 sample … 2 c 0.501 0.0199 Inf 0.462 0.540
## 3 random … 4 c 0.0399 0.0177 Inf 0.00523 0.0746
## 4 sample … 2 f 0.254 0.0171 Inf 0.220 0.287
## 5 sample … 2 f 0.311 0.0171 Inf 0.278 0.345
## 6 random … 10 f 0.0340 0.0171 Inf 0.000543 0.0675
## # … with 2 more variables: z.ratio <dbl>, p.value <dbl>Predictive value smallest pairwise effect commercial-lab
predictive.pairwise %>% dplyr::filter(p.value < .05, z.ratio > 0) %>% dplyr::group_by(grain, contrast) %>% filter(z.ratio == min(z.ratio))## # A tibble: 5 x 10
## # Groups: grain, contrast [5]
## contrast sampleSize grain estimate SE df asymp.LCL asymp.UCL
## <fct> <fct> <fct> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 sample … 2 c 0.737 0.0398 Inf 0.659 0.815
## 2 sample … 2 c 0.914 0.0422 Inf 0.831 0.997
## 3 random … 4 c 0.0978 0.0377 Inf 0.0240 0.172
## 4 sample … 2 f 0.439 0.0365 Inf 0.368 0.511
## 5 sample … 2 f 0.516 0.0365 Inf 0.444 0.587
## # … with 2 more variables: z.ratio <dbl>, p.value <dbl>Detection value smallest pairwise effect commercial-lab
detection.pairwise %>% dplyr::filter(p.value < .05, z.ratio > 0) %>% dplyr::group_by(grain, contrast) %>% filter(z.ratio == min(z.ratio))## # A tibble: 6 x 10
## # Groups: grain, contrast [6]
## contrast sampleSize grain estimate SE df asymp.LCL asymp.UCL
## <fct> <fct> <fct> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 sample … 2 c 0.343 0.0187 Inf 0.306 0.380
## 2 sample … 2 c 0.501 0.0199 Inf 0.462 0.540
## 3 random … 4 c 0.0399 0.0177 Inf 0.00523 0.0746
## 4 sample … 2 f 0.254 0.0171 Inf 0.220 0.287
## 5 sample … 2 f 0.311 0.0171 Inf 0.278 0.345
## 6 random … 10 f 0.0340 0.0171 Inf 0.000543 0.0675
## # … with 2 more variables: z.ratio <dbl>, p.value <dbl>Predictive value polynomial fit commercial-lab
predictive.poly[predictive.poly$testCond == "sample",]## contrast grain testCond estimate SE df asymp.LCL asymp.UCL
## 13 linear c sample 33.58800 0.9738438 Inf 31.67930 35.49670
## 14 quadratic c sample -26.37411 2.0076247 Inf -30.30898 -22.43924
## 15 cubic c sample 90.12431 26.5866237 Inf 38.01549 142.23314
## 16 linear f sample 27.82441 0.9482521 Inf 25.96587 29.68295
## 17 quadratic f sample -30.80938 1.9440620 Inf -34.61968 -26.99909
## 18 cubic f sample 237.24458 25.8205091 Inf 186.63731 287.85184
## z.ratio p.value
## 13 34.490135 1.127637e-260
## 14 -13.136972 2.021558e-39
## 15 3.389837 6.993427e-04
## 16 29.342840 2.948948e-189
## 17 -15.847943 1.452581e-56
## 18 9.188222 3.993872e-20Detection value polynomial fit commercial-lab
detection.poly[detection.poly$testCond == "sample",]## contrast grain testCond estimate SE df asymp.LCL asymp.UCL
## 13 linear c sample 16.80430 0.4577924 Inf 15.90704 17.70156
## 14 quadratic c sample -13.06385 0.9437606 Inf -14.91359 -11.21411
## 15 cubic c sample 43.90222 12.4980363 Inf 19.40652 68.39792
## 16 linear f sample 15.92080 0.4457602 Inf 15.04713 16.79448
## 17 quadratic f sample -17.42398 0.9138766 Inf -19.21514 -15.63281
## 18 cubic f sample 127.27518 12.1378629 Inf 103.48541 151.06495
## z.ratio p.value
## 13 36.707249 5.595930e-295
## 14 -13.842334 1.415362e-43
## 15 3.512729 4.435288e-04
## 16 35.716072 2.225900e-279
## 17 -19.066008 4.838610e-81
## 18 10.485798 1.003944e-25#Extra analyses
setwd('../data/raw')
ts.data.commercial <- read.delim("esMethods_Commercial_TimeSeries.txt", head = TRUE)
ts.data.evAct <- read.delim("esMethods_EvAct_TimeSeries.txt", head = TRUE)
names(ts.data.commercial)[names(ts.data.commercial) == "movName"] <- "movieName"
names(ts.data.evAct)[names(ts.data.evAct) == "movName"] <- "movieName"Correlate correlation between segmentation of different movies
comm.gp.bp <- ts.data.commercial %>% dplyr::group_by(grain, movieName, timeSeries) %>% dplyr::summarise(bp = mean(subject.bp))
evAct.gp.bp <- ts.data.evAct %>% dplyr::group_by(grain, movieName, timeSeries) %>% dplyr::summarise(bp = mean(subject.bp))
comm.corr.fine <- cor.test(comm.gp.bp$bp[comm.gp.bp$movieName == "3Iron" & comm.gp.bp$grain == "f"],
comm.gp.bp$bp[comm.gp.bp$movieName == "Corn" & comm.gp.bp$grain == "f"])
comm.corr.coarse <- cor.test(comm.gp.bp$bp[comm.gp.bp$movieName == "3Iron" & comm.gp.bp$grain == "c"],
comm.gp.bp$bp[comm.gp.bp$movieName == "Corn" & comm.gp.bp$grain == "c"])
evAct.corr.fine.bed_coffee <- cor.test(evAct.gp.bp$bp[evAct.gp.bp$movieName == "w_bed.mp4" & evAct.gp.bp$grain == "f"],
evAct.gp.bp$bp[evAct.gp.bp$movieName == "w_coffee.mp4" & evAct.gp.bp$grain == "f"])
evAct.corr.fine.bed_dishes <- cor.test(evAct.gp.bp$bp[evAct.gp.bp$movieName == "w_bed.mp4" & evAct.gp.bp$grain == "f"],
evAct.gp.bp$bp[evAct.gp.bp$movieName == "w_dishes.mp4" & evAct.gp.bp$grain == "f"])
evAct.corr.fine.bed_laundry <- cor.test(evAct.gp.bp$bp[evAct.gp.bp$movieName == "w_bed.mp4" & evAct.gp.bp$grain == "f"],
evAct.gp.bp$bp[evAct.gp.bp$movieName == "w_laundry.mp4" & evAct.gp.bp$grain == "f"])
evAct.corr.coarse.bed_coffee <- cor.test(evAct.gp.bp$bp[evAct.gp.bp$movieName == "w_bed.mp4" & evAct.gp.bp$grain == "c"],
evAct.gp.bp$bp[evAct.gp.bp$movieName == "w_coffee.mp4" & evAct.gp.bp$grain == "c"])
evAct.corr.coarse.bed_dishes <- cor.test(evAct.gp.bp$bp[evAct.gp.bp$movieName == "w_bed.mp4" & evAct.gp.bp$grain == "c"],
evAct.gp.bp$bp[evAct.gp.bp$movieName == "w_dishes.mp4" & evAct.gp.bp$grain == "c"])
evAct.corr.coarse.bed_laundry <- cor.test(evAct.gp.bp$bp[evAct.gp.bp$movieName == "w_bed.mp4" & evAct.gp.bp$grain == "c"],
evAct.gp.bp$bp[evAct.gp.bp$movieName == "w_laundry.mp4" & evAct.gp.bp$grain == "c"])
evAct.corr.fine.coffee_dishes <- cor.test(evAct.gp.bp$bp[evAct.gp.bp$movieName == "w_coffee.mp4" & evAct.gp.bp$grain == "f"],
evAct.gp.bp$bp[evAct.gp.bp$movieName == "w_dishes.mp4" & evAct.gp.bp$grain == "f"])
evAct.corr.fine.coffee_laundry <- cor.test(evAct.gp.bp$bp[evAct.gp.bp$movieName == "w_coffee.mp4" & evAct.gp.bp$grain == "f"],
evAct.gp.bp$bp[evAct.gp.bp$movieName == "w_laundry.mp4" & evAct.gp.bp$grain == "f"])
evAct.corr.coarse.coffee_dishes <- cor.test(evAct.gp.bp$bp[evAct.gp.bp$movieName == "w_coffee.mp4" & evAct.gp.bp$grain == "c"],
evAct.gp.bp$bp[evAct.gp.bp$movieName == "w_dishes.mp4" & evAct.gp.bp$grain == "c"])
evAct.corr.coarse.coffee_laundry <- cor.test(evAct.gp.bp$bp[evAct.gp.bp$movieName == "w_coffee.mp4" & evAct.gp.bp$grain == "c"],
evAct.gp.bp$bp[evAct.gp.bp$movieName == "w_laundry.mp4" & evAct.gp.bp$grain == "c"])
evAct.corr.fine.dishes_laundry <- cor.test(evAct.gp.bp$bp[evAct.gp.bp$movieName == "w_dishes.mp4" & evAct.gp.bp$grain == "f"],
evAct.gp.bp$bp[evAct.gp.bp$movieName == "w_laundry.mp4" & evAct.gp.bp$grain == "f"])
evAct.corr.coarse.dishes_laundry <- cor.test(evAct.gp.bp$bp[evAct.gp.bp$movieName == "w_dishes.mp4" & evAct.gp.bp$grain == "c"],
evAct.gp.bp$bp[evAct.gp.bp$movieName == "w_laundry.mp4" & evAct.gp.bp$grain == "c"])